source: trunk/MagicSoft/Mars/mpedestal/MPedestalSubtractedEvt.cc@ 8339

/* ======================================================================== *\
! $Name: not supported by cvs2svn $:$Id: MPedestalSubtractedEvt.cc,v 1.2 2007-01-15 12:06:15 tbretz Exp $
! --------------------------------------------------------------------------
!
! *
! * 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): Thomas Bretz, 10/2006 <mailto:tbretz@astro.uni-wuerzburg.de>
!
!   Copyright: MAGIC Software Development, 2000-2006
!
!
\* ======================================================================== */

/////////////////////////////////////////////////////////////////////////////
//
//  MPedestalSubtractedEvt
//
//  Storage container to store the raw FADC values.
//
/////////////////////////////////////////////////////////////////////////////

#include "MPedestalSubtractedEvt.h"

ClassImp(MPedestalSubtractedEvt);

using namespace std;

// --------------------------------------------------------------------------
//
// Initialize number of samples (per pixel) and number of pixels.
//
// Initialize the correct length of fSamples and fSamplesRaw
//
// And reset its contents to 0.
//
void MPedestalSubtractedEvt::InitSamples(UInt_t samples, UInt_t pixels)
{
    fNumSamples = samples;

    if (pixels>0)
        fNumPixels = pixels;

    fSamples.Set(fNumPixels*fNumSamples);
    fSamplesRaw.Set(fNumPixels*fNumSamples);

    fSamples.Reset();
    fSamplesRaw.Reset();
}

// --------------------------------------------------------------------------
//
// Return a pointer to the first slice with subtracted pedestal of
// the samples of the pixel with given pixel-index. If the number
// exceeds the number of pixels NULL is returned.
//
// The user is responsible not to exceed the slices for one pixel!
//
Float_t *MPedestalSubtractedEvt::GetSamples(UInt_t pixel) const
{
    return pixel>=fNumPixels ? NULL : fSamples.GetArray()+pixel*fNumSamples;
}

// --------------------------------------------------------------------------
//
// Return a pointer to the first slice of the raw-data samples of the pixel
// with given pixel-index. If the number exceeds the number of
// pixels NULL is returned.
//
// The user is responsible not to exceed the slices for one pixel!
//
Byte_t *MPedestalSubtractedEvt::GetSamplesRaw(UInt_t pixel) const
{
    return pixel>=fNumPixels ? NULL : fSamplesRaw.GetArray()+pixel*fNumSamples;
}

// --------------------------------------------------------------------------
//
// Return some information about saturation in the raw-data of pixel idx.
//
// The search range is defined by [first,last]. Saturation is considered if
// contents is >= limit.
//
// The number of saturating slices are returned and first/last are filled
// with the first and last saturating slice index w.r.t. the beginning of
// the raw-data of this pixel not first.
//
// Warning: No range checks and no sanity checks are done!
//
Int_t MPedestalSubtractedEvt::GetSaturation(const Int_t idx, Int_t limit, Int_t &first, Int_t &last) const
{
    // Determin saturation of hi-gains
    Byte_t *p0 = GetSamplesRaw(idx);

    Byte_t *sathi0 = 0; // first saturating hi-gain slice
    Byte_t *sathi1 = 0; // last  saturating hi-gain slice

    Int_t num = 0;

    const Byte_t *end = p0+last;
    for (Byte_t *ptr=p0+first; ptr<=end; ptr++)
    {
        if (*ptr>=limit)
        {
            sathi1 = ptr;
            if (!sathi0)
                sathi0 = ptr;
            num++;
        }
    }

    last  = sathi1 ? sathi1-p0 : -1;
    first = sathi0 ? sathi0-p0 : -1;

    return num;
}

// --------------------------------------------------------------------------
//
// Get the maximum of the slices [first,last] of pixel index.
//
// The position returned is the index of the position of the pedestal
// subtracted maximum w.r.t. to first.
// The value returned is the maximum of the raw-data.
//
// Warning: No range checks and no sanity checks are done!
//
Int_t MPedestalSubtractedEvt::GetMax(const Int_t idx, const Int_t first, const Int_t last, Int_t &val) const
{
    // Get pointer to first slice to be considered
    Byte_t  const *samb = GetSamplesRaw(idx);
    Float_t const *samf = GetSamples(idx);

    Byte_t  const *ptrb = samb+first;
    Float_t const *ptrf = samf+first;

    // The best information so far: the first slice is the maximum
    const Byte_t  *maxb = ptrb;
    const Float_t *maxf = ptrf;

    // Store the pointer to the first slice
//    const Byte_t  *begb = ptrb;
    const Float_t *begf = ptrf;

    // Calculate the last slice to be considered
    const Byte_t  *endb = samb+last;

    while (ptrb<endb)
    {
        // Pre-increment: check the second slice first
        if (*++ptrb>*maxb)
            maxb = ptrb;
        if (*++ptrf>*maxf)
            maxf = ptrf;
    }

    val = *maxb;
    return maxf-begf;
}

/*
#include <TSpline.h>
#include "MArrayD.h"
void  MPedestalSubtractedEvt::InterpolateSaturation(const Int_t idx, Int_t limit, Int_t first, Int_t last) const
{
    MArrayD x(GetNumSamples());
    MArrayD y(GetNumSamples());

    Float_t *s = GetSamples(idx);

    Int_t n = 0;
    for (unsigned int i=0; i<GetNumSamples(); i++)
    {
        if (s[i]>limit)
            continue;
        x[n] = i;
        y[n] = s[i];
        n++;
    }

    TSpline5 sp("", x.GetArray(), y.GetArray(), n);

    for (unsigned int i=0; i<GetNumSamples(); i++)
    {
        if (s[i]>limit)
            s[i] = sp.Eval(i);
    }
}
*/