Index: trunk/MagicSoft/Mars/msignal/MExtractFixedWindowPeakSearch.cc
===================================================================
--- trunk/MagicSoft/Mars/msignal/MExtractFixedWindowPeakSearch.cc	(revision 3862)
+++ trunk/MagicSoft/Mars/msignal/MExtractFixedWindowPeakSearch.cc	(revision 3862)
@@ -0,0 +1,332 @@
+/* ======================================================================== *\
+   !
+   ! *
+   ! * This file is part of MARS, the MAGIC Analysis and Reconstruction
+   ! * Software. It is distributed to you in the hope that it can be a useful
+   ! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
+   ! * It is distributed WITHOUT ANY WARRANTY.
+   ! *
+   ! * Permission to use, copy, modify and distribute this software and its
+   ! * documentation for any purpose is hereby granted without fee,
+   ! * provided that the above copyright notice appear in all copies and
+   ! * that both that copyright notice and this permission notice appear
+   ! * in supporting documentation. It is provided "as is" without express
+   ! * or implied warranty.
+   ! *
+   !
+   !
+   !   Author(s): Abelardo Moralejo, 4/2004 <mailto:moralejo@pd.infn.it>
+   !              Markus Gaug      , 4/2004 <mailto:markus@ifae.es>
+   !   Copyright: MAGIC Software Development, 2000-2004
+   !
+   !
+   \* ======================================================================== */
+
+//////////////////////////////////////////////////////////////////////////////
+//
+//   MExtractFixedWindowPeakSearch
+//
+//  Calculate the signal integrating fWindowSize time slices, the same for 
+//  all pixels. The integrated slices change from event to event, since the
+//  pulse positions in the FADC jump between events, but apparently in a 
+//  "coherent" fashion. We first loop over all pixels and find the one 
+//  which has the highest sum of fPeakSearchWindowSize (default: 4) consecutive 
+//  non-saturated high gain slices. The position of the peak in this pixel 
+//  determines the integration window position for all pixels of this event. 
+//  For the moment we neglect time delays between pixels (which are in most 
+//  cases small). The class is based on MExtractSignal2.
+//
+//////////////////////////////////////////////////////////////////////////////
+#include "MExtractFixedWindowPeakSearch.h"
+#include "MExtractor.h"
+
+#include "MLog.h"
+#include "MLogManip.h"
+
+#include "MRawEvtPixelIter.h"
+
+#include "MExtractedSignalCam.h"
+#include "MExtractedSignalPix.h"
+
+#include "MPedestalCam.h"
+#include "MPedestalPix.h"
+
+ClassImp(MExtractFixedWindowPeakSearch);
+
+using namespace std;
+
+const Byte_t MExtractFixedWindowPeakSearch::fgHiGainWindowSize     = 6;
+const Byte_t MExtractFixedWindowPeakSearch::fgLoGainWindowSize     = 6;
+const Byte_t MExtractFixedWindowPeakSearch::fgPeakSearchWindowSize = 4;
+const Byte_t MExtractFixedWindowPeakSearch::fgHiGainFirst          = 0;
+const Byte_t MExtractFixedWindowPeakSearch::fgHiGainLast           = 14;
+const Byte_t MExtractFixedWindowPeakSearch::fgLoGainFirst          = 0;
+const Byte_t MExtractFixedWindowPeakSearch::fgLoGainLast           = 14;
+// --------------------------------------------------------------------------
+//
+// Default constructor. 
+//
+MExtractFixedWindowPeakSearch::MExtractFixedWindowPeakSearch(const char *name, const char *title)
+{
+
+  fName  = name  ? name  : "MExtractFixedWindowPeakSearch";
+  fTitle = title ? title : "Task to extract the signal from the FADC slices";
+
+  SetRange(fgHiGainFirst, fgHiGainLast, fgLoGainFirst, fgLoGainLast);
+  SetWindows();
+}
+
+void MExtractFixedWindowPeakSearch::SetWindows(Byte_t windowh, Byte_t windowl, Byte_t peaksearchwindow)
+{
+  //
+  // Set windows to even number of slices due to clock noise (odd/even slice effect).
+  //
+  fWindowSizeHiGain = windowh & ~1;
+  fWindowSizeLoGain = windowl & ~1;
+  fPeakSearchWindowSize = peaksearchwindow & ~1;
+
+
+  if (fWindowSizeHiGain != windowh)
+    *fLog << endl << warn << 
+      "MExtractFixedWindowPeakSearch::SetWindows - Hi Gain window size has to be even, set to: " 
+	  << int(fWindowSizeHiGain) << " samples " << endl;
+    
+  if (fWindowSizeLoGain != windowl)
+    *fLog << endl << warn << 
+      "MExtractFixedWindowPeakSearch::SetWindows - Lo Gain window size has to be even, set to: " 
+	  << int(fWindowSizeLoGain) << " samples " << endl;
+
+  if (fPeakSearchWindowSize != peaksearchwindow)
+    *fLog << endl << warn << 
+      "MExtractFixedWindowPeakSearch::SetWindows - Peak Search window size has to be even, set to: " 
+	  << int(fPeakSearchWindowSize) << " samples " << endl;
+
+
+  if (fWindowSizeHiGain<2) 
+    {
+      fWindowSizeHiGain = 2;
+      *fLog << warn 
+            << "MExtractFixedWindowPeakSearch::SetWindows - Hi Gain window size set to two samples" << endl;
+    }
+
+  if (fWindowSizeLoGain<2) 
+    {
+      fWindowSizeLoGain = 2;
+      *fLog << warn 
+            << "MExtractFixedWindowPeakSearch::SetWindows - Lo Gain window size set to two samples" << endl;
+    }
+
+  if (fPeakSearchWindowSize<2) 
+    {
+      fPeakSearchWindowSize = 2;
+      *fLog << warn 
+            << "MExtractFixedWindowPeakSearch::SetWindows - Peak Search window size set to two samples" << endl;
+    }
+
+
+  fNumHiGainSamples = (Float_t)fWindowSizeHiGain;
+  fNumLoGainSamples = (Float_t)fWindowSizeLoGain;
+
+  fSqrtHiGainSamples = TMath::Sqrt(fNumHiGainSamples);
+  fSqrtLoGainSamples = TMath::Sqrt(fNumLoGainSamples);
+}
+
+// --------------------------------------------------------------------------
+//
+// FindPeak
+// Finds highest sum of "window" consecutive FADC slices in a pixel, and store
+// in "startslice" the first slice of the group which yields the maximum sum.
+// In "max" the value of the maximum sum is stored, in "sat" the number of 
+// saturated slices.
+//
+void MExtractFixedWindowPeakSearch::FindPeak(Byte_t *ptr, Byte_t window, Byte_t &startslice, 
+			       Int_t &max, Int_t &sat) const
+{
+  const Byte_t *end = ptr + fHiGainLast - fHiGainFirst + 1;
+
+  Int_t sum=0;
+
+  //
+  // Calculate the sum of the first "window" slices
+  //
+  sat = 0;
+  Byte_t *p = ptr;
+
+  while (p<ptr+window)
+    {
+      sum += *p;
+      if (*p++ >= fSaturationLimit)
+	sat++;
+    }
+
+  //
+  // Check for saturation in all other slices
+  //
+  while (p<end)
+    if (*p++ >= fSaturationLimit)
+      sat++;
+
+  //
+  // Calculate the i-th sum of n as
+  //    sum_i+1 = sum_i + slice[i+window] - slice[i]
+  // This is fast and accurate (because we are using int's)
+  //
+  max=sum;
+  for (p=ptr; p+window<end; p++)
+    {
+      sum += *(p+window) - *p;
+      if (sum > max)
+	{
+	  max = sum;
+	  startslice = p-ptr;
+	}
+    }
+
+  return;
+}
+
+// --------------------------------------------------------------------------
+//
+// FindSignal
+// Adds up "window" slices starting in slice to which "ptr" points. The result 
+// is stored in "sum", and the number of saturated  slices in "sat".
+//
+void MExtractFixedWindowPeakSearch::FindSignalHiGain(Byte_t *ptr, Int_t size, Int_t &sum, Byte_t &sat) const
+{
+  //
+  // Calculate the sum of the "window" slices starting in ptr
+  //
+  sum = 0;
+  sat = 0;
+  Byte_t *p = ptr;
+
+  while (p<ptr+fWindowSizeHiGain)
+    {
+      sum += *p;
+      if (*p++ >= fSaturationLimit)
+	sat++;
+    }
+
+  return;
+}
+
+// --------------------------------------------------------------------------
+//
+// FindSignal
+// Adds up "window" slices starting in slice to which "ptr" points. The result 
+// is stored in "sum", and the number of saturated  slices in "sat".
+//
+void MExtractFixedWindowPeakSearch::FindSignalLoGain(Byte_t *ptr, Int_t size, Int_t &sum, Byte_t &sat) const
+{
+  //
+  // Calculate the sum of the "window" slices starting in ptr
+  //
+  sum = 0;
+  sat = 0;
+  Byte_t *p = ptr;
+
+  while (p<ptr+fWindowSizeLoGain)
+    {
+      sum += *p;
+      if (*p++ >= fSaturationLimit)
+	sat++;
+    }
+
+  return;
+}
+
+// --------------------------------------------------------------------------
+//
+// Process
+// First find the pixel with highest sum of fPeakSearchWindowSize slices (default:4)
+// "startslice" will mark the slice at which the highest sum begins for that pixel.
+// Then define the beginning of the integration window for ALL pixels as the slice
+// before that: startslice-1 (this is somehow arbitrary), unless of course startslice=0,
+// in which case we start at 0. We will also check that the integration window does not 
+// go beyond the FADC limits.
+//
+Int_t MExtractFixedWindowPeakSearch::Process()
+{
+
+  MRawEvtPixelIter pixel(fRawEvt);
+
+  Int_t sat;
+  Int_t maxsumhi = 0;
+  Byte_t startslice;
+  Byte_t hiGainFirst = 0;
+  Byte_t loGainFirst = 0;
+
+  while (pixel.Next())
+    {
+      Int_t sumhi;
+      sat = 0;
+
+      FindPeak(pixel.GetHiGainSamples(), fPeakSearchWindowSize, startslice, sumhi, sat);
+      if (sumhi > maxsumhi && sat == 0 )
+	{
+	  maxsumhi = sumhi;
+	  if (startslice > 0)
+	    hiGainFirst = startslice-1;
+	  else
+	    hiGainFirst = 0;
+	}
+    }
+
+
+  loGainFirst = hiGainFirst;
+
+  // Make sure we will not integrate beyond the hi gain limit:
+  if (hiGainFirst+fWindowSizeHiGain > pixel.GetNumHiGainSamples())
+    hiGainFirst = pixel.GetNumHiGainSamples()-fWindowSizeHiGain;
+
+  // Make sure we will not integrate beyond the lo gain limit:
+  if (loGainFirst+fWindowSizeLoGain > pixel.GetNumLoGainSamples())
+    loGainFirst = pixel.GetNumLoGainSamples()-fWindowSizeLoGain;
+
+  pixel.Reset();
+  fSignals->Clear();
+
+  sat = 0;
+  while (pixel.Next())
+    {
+      //
+      // Find signal in hi- and lo-gain
+      //
+      Int_t sumhi;
+      Byte_t sathi;
+
+      FindSignalHiGain(pixel.GetHiGainSamples()+hiGainFirst, fWindowSizeHiGain, sumhi, sathi);
+
+      Int_t sumlo=0;
+      Byte_t satlo=0;
+      if (pixel.HasLoGain())
+        {
+	  FindSignalLoGain(pixel.GetLoGainSamples()+loGainFirst, fWindowSizeLoGain, sumlo, satlo);
+	  if (satlo)
+	    sat++;
+        }
+
+      //
+      // Take corresponding pedestal
+      //
+      const Int_t pixid = pixel.GetPixelId();
+
+      const MPedestalPix  &ped = (*fPedestals)[pixid];
+      MExtractedSignalPix &pix = (*fSignals)[pixid];
+
+      const Float_t pedes  = ped.GetPedestal();
+      const Float_t pedrms = ped.GetPedestalRms();
+
+      //
+      // Set extracted signal with pedestal substracted
+      //
+      pix.SetExtractedSignal(sumhi - pedes*fNumHiGainSamples, pedrms*fSqrtHiGainSamples,
+                             sumlo - pedes*fNumLoGainSamples, pedrms*fSqrtLoGainSamples);
+
+      pix.SetGainSaturation(sathi, sathi, satlo);
+    } /* while (pixel.Next()) */
+
+  fSignals->SetReadyToSave();
+
+  return kTRUE;
+}
Index: trunk/MagicSoft/Mars/msignal/Makefile
===================================================================
--- trunk/MagicSoft/Mars/msignal/Makefile	(revision 3861)
+++ trunk/MagicSoft/Mars/msignal/Makefile	(revision 3862)
@@ -37,4 +37,8 @@
 SRCFILES = MExtractedSignalCam.cc \
            MExtractedSignalPix.cc \
+           MExtractor.cc \
+           MExtractFixedWindow.cc \
+           MExtractSlidingWindow.cc \
+           MExtractFixedWindowPeakSearch.cc \
            MExtractSignal.cc \
            MExtractSignal2.cc \
Index: trunk/MagicSoft/Mars/msignal/SignalLinkDef.h
===================================================================
--- trunk/MagicSoft/Mars/msignal/SignalLinkDef.h	(revision 3861)
+++ trunk/MagicSoft/Mars/msignal/SignalLinkDef.h	(revision 3862)
@@ -10,4 +10,10 @@
 #pragma link C++ class MExtractSignal2+;
 #pragma link C++ class MExtractSignal3+;
+
+#pragma link C++ class MExtractor+;
+#pragma link C++ class MExtractFixedWindow+;
+#pragma link C++ class MExtractSlidingWindow+;
+#pragma link C++ class MExtractFixedWindowPeakSearch+;
+   
 
 #pragma link C++ class MExtractPINDiode++;
