Index: /trunk/FACT++/src/DrsCalib.h
===================================================================
--- /trunk/FACT++/src/DrsCalib.h	(revision 11678)
+++ /trunk/FACT++/src/DrsCalib.h	(revision 11678)
@@ -0,0 +1,261 @@
+class CalibData
+{
+protected:
+    uint64_t fNumEntries;
+    uint64_t fNumCounter;
+
+    size_t fNumSamples;
+    size_t fNumChannels;
+
+    vector<int64_t> fSum;
+    vector<int64_t> fSum2;
+
+    /*
+    vector<map<int16_t, uint8_t> > fMedian;
+    vector<int16_t> fResult;
+     */
+public:
+    CalibData() : fNumEntries(0), fNumSamples(0), fNumChannels(0) { }
+    void Reset()
+    {
+        fNumEntries  = 0;
+        fNumSamples  = 0;
+        fNumChannels = 0;
+
+        fSum.clear();
+        fSum2.clear();
+    }
+
+    void InitSize(uint16_t channels, uint16_t samples)
+    {
+        fNumChannels = channels;
+        fNumSamples  = samples;
+
+        fSum.resize(samples*channels);
+        fSum2.resize(samples*channels);
+
+        //fMedian.resize(40*9*1024);
+    }
+/*
+                continue;
+
+                if (ch<40*9)
+                {
+                    fMedian[abs][val[rel]]++;
+
+                    int n= 8;
+                    if (fNumEntries>0 && fNumEntries%(n*100)==0)
+                    {
+                        fResult.resize(40*9*1024);
+
+                        uint16_t sum = 0;
+                        for (map<int16_t, uint8_t>::const_iterator it=fMedian[abs].begin();
+                             it!=fMedian[abs].end(); it++)
+                        {
+                            map<int16_t, uint8_t>::const_iterator is = it;
+                            is++;
+                            if (is==fMedian[abs].end())
+                                break;
+
+                            sum += it->second;
+
+                            double med = 0;
+
+                            med += int64_t(it->second)*int64_t(it->first);
+                            med += int64_t(is->second)*int64_t(is->first);
+                            med /= int64_t(it->second)+int64_t(is->second);
+
+                            if (sum==n*50)
+                            {
+                                fResult[abs] = it->first;
+                                cout << ch << " MED+(50): " << it->first << endl;
+                            }
+                            if (sum<n*50 && sum+is->second>n*50)
+                            {
+                                fResult[abs] = med;
+                                cout << ch << " MED-(50): " << med << endl;
+                            }
+                        }
+                        cout << ch << " AVG=" << float(fSum[abs])/fNumEntries << endl;
+                        fMedian[abs].clear();
+                    }
+                } // -2029 -2012 -2003 -1996 -1990 // -1955
+                */
+
+    void AddRel(const int16_t *val, const int16_t *start)
+    {
+        for (size_t ch=0; ch<fNumChannels; ch++)
+        {
+            const int16_t spos = start[ch];
+            if (spos<0)
+                continue;
+
+            const size_t pos = ch*fNumSamples;
+
+            for (size_t i=0; i<fNumSamples; i++)
+            {
+                // Value is relative to trigger
+                // Abs is corresponding index relative to DRS pipeline
+                const size_t rel = pos +  i;
+                const size_t abs = pos + (spos+i)%fNumSamples;
+
+                const int64_t v = val[rel];
+
+                fSum[abs]  += v;
+                fSum2[abs] += v*v;
+            }
+        }
+
+        fNumEntries++;
+    }
+
+    void AddRel(const int16_t *val,    const int16_t *start,
+                const int32_t *offset, const uint32_t scale)
+    {
+        for (size_t ch=0; ch<fNumChannels; ch++)
+        {
+            const int16_t spos = start[ch];
+            if (spos<0)
+                continue;
+
+            const size_t pos = ch*fNumSamples;
+
+            for (size_t i=0; i<fNumSamples; i++)
+            {
+                // Value is relative to trigger
+                // Offset is relative to DRS pipeline
+                // Abs is corresponding index relative to DRS pipeline
+                const size_t rel = pos +  i;
+                const size_t abs = pos + (spos+i)%fNumSamples;
+
+                const int64_t v = int64_t(val[rel])*scale-offset[abs];
+
+                fSum[abs]  += v;
+                fSum2[abs] += v*v;
+            }
+        }
+
+        fNumEntries++;
+    }
+
+    void AddAbs(const int16_t *val,    const  int16_t *start,
+                const int32_t *offset, const uint32_t scale)
+    {
+        for (size_t ch=0; ch<fNumChannels; ch++)
+        {
+            const int16_t spos = start[ch];
+            if (spos<0)
+                continue;
+
+            const size_t pos = ch*fNumSamples;
+
+            for (size_t i=0; i<fNumSamples; i++)
+            {
+                // Value is relative to trigger
+                // Offset is relative to DRS pipeline
+                // Abs is corresponding index relative to DRS pipeline
+                const size_t rel = pos +  i;
+                const size_t abs = pos + (spos+i)%fNumSamples;
+
+                const int64_t v = int64_t(val[rel])*scale-offset[abs];
+
+                fSum[rel]  += v;
+                fSum2[rel] += v*v;
+            }
+        }
+
+        fNumEntries++;
+    }
+
+    static void Apply(int16_t *val, const int16_t *start, uint32_t roi,
+                      const int32_t *offset, const uint32_t scaleabs,
+                      const int32_t *gain,   const uint32_t scalegain,
+                      const int32_t *trgoff, const uint32_t scalerel)
+    {
+        for (size_t ch=0; ch<1440; ch++)
+        {
+            const int16_t spos = start[ch];
+            if (spos<0)
+                continue;
+
+            const size_t pos = ch*roi;
+
+            for (size_t i=0; i<roi; i++)
+            {
+                // Value is relative to trigger
+                // Offset is relative to DRS pipeline
+                // Abs is corresponding index relative to DRS pipeline
+                const size_t rel = pos +  i;
+                const size_t abs = pos + (spos+i)%1024;
+
+                const int64_t v =
+                    + int64_t(val[rel])   *scaleabs*scalerel
+                    - int64_t(trgoff[rel])*scaleabs
+                    - int64_t(offset[abs])*scalerel
+                    ;
+
+                const int64_t div = int64_t(gain[abs])*scaleabs*scalerel;
+
+                val[rel] = (v*scalegain)/div;
+            }
+        }
+    }
+
+    pair<vector<double>,vector<double> > GetSampleStats() const
+    {
+        if (fNumEntries==0)
+            return make_pair(vector<double>(),vector<double>());
+
+        vector<double> mean(fSum.size());
+        vector<double> error(fSum.size());
+
+        vector<int64_t>::const_iterator it = fSum.begin();
+        vector<int64_t>::const_iterator i2 = fSum2.begin();
+        vector<double>::iterator        im = mean.begin();
+        vector<double>::iterator        ie = error.begin();
+
+        int cnt = 0;
+        while (it!=fSum.end())
+        {
+            *im = /*cnt<fResult.size() ? fResult[cnt] :*/ double(*it)/fNumEntries;
+            *ie = sqrt(double(*i2*int64_t(fNumEntries) - *it * *it))/fNumEntries;
+
+            im++;
+            ie++;
+            it++;
+            i2++;
+            cnt++;
+        }
+
+
+        /*
+         valarray<double> ...
+
+         mean /= fNumEntries;
+         error = sqrt(error/fNumEntries - mean*mean);
+         */
+
+        return make_pair(mean, error);
+    }
+
+    const vector<int64_t> &GetSum() const { return fSum; }
+
+    uint64_t GetNumEntries() const { return fNumEntries; }
+};
+
+
+/*
+Int_t TH1::MyFill(Int_t x, Double_t w)
+{
+   fEntries++;
+
+   fArrayBin[x] += w;
+
+   fSumw2.fArray[x] += w*w
+   fTsumw   += w;
+   fTsumw2  += w*w;
+   fTsumwx  += w*x;
+   fTsumwx2 += w*x*x;
+   return bin;
+}
+*/