Ignore:
Timestamp:
01/06/04 16:12:03 (21 years ago)
Author:
raducci
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars/manalysis
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • trunk/MagicSoft/Mars/manalysis/MArrivalTime.cc

    r2638 r2752  
    2727// MArrivalTime
    2828//
    29 // P R E L I M I N A R Y
    30 // Do not use this container. It has yet to be defined.
     29// Times are calculated using the TSpline5 Root Class
     30//
    3131/////////////////////////////////////////////////////////////////////////////
    3232#include "MArrivalTime.h"
     
    5555
    5656//
    57 // Calculates the arrival time for each pixel (for now, simply by finding the peak)
     57// Calculates the arrival time for each pixel
     58// Possible Methods
     59// Case 1: Spline5 (From TSpline5 Root Class)
     60//
    5861//
    5962
     
    6669
    6770    MRawEvtPixelIter pixel((MRawEvtData*)&evt);
    68 
    69     Int_t saturatedpixels = 0;
    70 
     71   
    7172    while ( pixel.Next() )
    7273    {
    73         const UInt_t idx = pixel.GetPixelId();
     74        const UInt_t idx = pixel.GetPixelId();
    7475
    75         Byte_t *ptr = pixel.GetHiGainSamples();
     76//If pixel has saturated and hasn't lo gains we return -1
    7677
    77         Int_t n = evt.GetNumHiGainSamples();
     78        if (pixel.GetMaxHiGainSample() == 0xff)
     79        {   
     80            if (!pixel.HasLoGain())       
     81            {
     82                fData[idx]=-1.0;
     83            }
     84            else //Calculate time using Lo Gains
     85            {
     86                Byte_t *ptr = pixel.GetLoGainSamples();
    7887
    79         Int_t i;
     88                Int_t nSlice = evt.GetNumLoGainSamples();
     89//Some casts are needed because TSpline5 constructor accepts only Double_t values
     90                Double_t ptr2[nSlice];
     91                for (Int_t i = 0; i < nSlice; i++)
     92                    ptr2[i]=(Double_t)ptr[i];
     93                TSpline5 *spline = new TSpline5("spline",(Double_t) 0,(Double_t)(nSlice - 1),ptr2,nSlice);
     94                Double_t abscissa=0.0;
     95                Double_t maxAb=0.0;
     96                Double_t maxOrd=0.0;
     97                Double_t swap;
     98                while (abscissa <= nSlice - 1)
     99                {
     100                    swap=spline->Eval(abscissa);
     101                    if (swap > maxOrd)
     102                    {
     103                        maxOrd = swap;
     104                        maxAb  = abscissa;
     105                    }
     106                    abscissa += 0.1;
     107                }
     108                fData[idx]=maxAb;
     109            }
     110        }
     111        else //Calculate time using Hi Gains
     112        {
     113            Byte_t *ptr = pixel.GetHiGainSamples();
    80114
    81         Double_t arrtime = 0;
    82 
    83         Double_t maxsign = 0;
    84 
    85         for (i=0; i<n; i++)
    86         {
    87             if (ptr[i]==0xff)
    88                 break;
    89             if (ptr[i]>maxsign)
    90                {
    91                 maxsign = ptr[i];
    92                 arrtime = i+1;
    93                }
     115            Int_t nSlice = evt.GetNumHiGainSamples();
     116//Some casts are needed because TSpline5 constructor accepts only Double_t values
     117            Double_t  ptr2[nSlice];
     118            for (Int_t i = 0; i < nSlice; i++)
     119                ptr2[i]=(Double_t)ptr[i];
     120            TSpline5 *spline = new TSpline5("spline",(Double_t) 0,(Double_t)(nSlice - 1),ptr2,nSlice);
     121            Double_t abscissa=0.0;
     122            Double_t maxAb=0.0;
     123            Double_t maxOrd=0.0;
     124            Double_t swap;
     125            while (abscissa <= nSlice - 1)
     126            {
     127                swap=spline->Eval(abscissa);
     128                if (swap > maxOrd)
     129                {
     130                    maxOrd = swap;
     131                    maxAb  = abscissa;
     132                }
     133                abscissa += 0.1;
     134            }
     135            fData[idx]=maxAb;
    94136        }
    95 
    96         Bool_t saturatedlg = kFALSE;
    97 
    98         if (i!=n)
    99         {
    100             arrtime = 0;
    101             maxsign = 0;
    102 
    103             ptr = pixel.GetLoGainSamples();
    104             if (ptr==NULL)
    105             {
    106              *fLog << warn << "WARNING - Pixel #" << idx
    107                    << " saturated but has no low gains... skipping!" << endl;
    108              return;
    109             }
    110 
    111             for (i=0; i<n; i++)
    112             {
    113                 if (ptr[i]==0xff)
    114                     saturatedlg = kTRUE;
    115 
    116                 if (ptr[i]>maxsign)
    117                    {
    118                     maxsign = ptr[i];
    119                     arrtime = i+1;
    120                    }
    121             }
    122         }
    123 
    124         fData[idx]=arrtime;
    125 
    126         if (saturatedlg)
    127           saturatedpixels++;
    128137    }
    129 
    130     if (saturatedpixels>0)
    131         *fLog << warn << "WARNING: " << saturatedpixels
    132               << " pixel(s) had saturating low gains..." << endl;
    133 
    134 
    135138}
    136139
     
    138141// --------------------------------------------------------------------------
    139142//
    140 // Returns the arrival time value (for now, it's the FADC slice number).
     143// Returns the arrival time value
    141144//
    142145
  • trunk/MagicSoft/Mars/manalysis/MArrivalTime.h

    r2651 r2752  
    22#define MARS_MArrivalTime
    33
    4 #ifndef ROOT_TArrayD
    5 #include <TArrayD.h>
     4#ifndef ROOT_TArrayF
     5#include <TArrayF.h>
     6#endif
     7#ifndef ROOT_TSpline
     8#include <TSpline.h>
    69#endif
    710#ifndef MARS_MCamEvent
     
    1619{
    1720private:
    18     TArrayD fData;  // Stores the arrival times
     21    TArrayF fData;  // Stores the arrival times
    1922
    2023public:
     
    2629    void Calc(const MRawEvtData &evt, const MGeomCam &geom); // Calculates arrival times
    2730
    28     const TArrayD &GetData() const { return fData; }
     31    const TArrayF &GetData() const { return fData; }
    2932
    3033    Double_t operator[](int i) { return fData[i]; }
  • trunk/MagicSoft/Mars/manalysis/MArrivalTimeCalc.cc

    r2659 r2752  
    2828//
    2929//   This is a task that calculates the arrival times of photons.
    30 //   For now, it returns the number of time slice containig the maximum value.
    31 //
    32 //   P R E L I M I N A R Y
    33 //   Other more sophisticated methods have to be implemented.
     30//   It returns the absolute maximum of the spline that interpolates
     31//   the FADC slices
    3432//
    3533// Input Containers:
     
    3735//
    3836// Output Containers:
    39 //   //MArrivalTime
     37//   MArrivalTime
    4038//   MRawEvtData
    4139//////////////////////////////////////////////////////////////////////////////
    4240
    43 //#include "MArrivalTime.h"
     41#include "MArrivalTime.h"
    4442#include "MArrivalTimeCalc.h"
    4543
     
    6563// Default constructor.
    6664//
     65
    6766MArrivalTimeCalc::MArrivalTimeCalc(const char *name, const char *title)
    6867{
     
    7069    fTitle = title ? title : "Calculate photons arrival time";
    7170
    72     AddToBranchList("MRawEvtData.fHiGainPixId");
    73     AddToBranchList("MRawEvtData.fLoGainPixId");
    74     AddToBranchList("MRawEvtData.fHiGainFadcSamples");
    75     AddToBranchList("MRawEvtData.fLoGainFadcSamples");
     71//    AddToBranchList("MRawEvtData.fHiGainPixId");
     72//    AddToBranchList("MRawEvtData.fLoGainPixId");
     73//    AddToBranchList("MRawEvtData.fHiGainFadcSamples");
     74//    AddToBranchList("MRawEvtData.fLoGainFadcSamples");
    7675
    7776}
     
    8281//  - MRawRunHeader
    8382//  - MRawEvtData
    84 //  //- MArrivalTime
     83//  - MArrivalTime
    8584//  - MGeomCam
    8685//
    8786// The following output containers are also searched and created if
    8887// they were not found:
    89 //  //- MArrivalTime
     88//  - MArrivalTime
    9089//  - MRawEvtData
    9190//
     
    114113    }
    115114
    116 /*    fArrTime = (MArrivalTime*)pList->FindCreateObj(AddSerialNumber("MArrivalTime"));
     115    fArrTime = (MArrivalTime*)pList->FindCreateObj(AddSerialNumber("MArrivalTime"));
    117116    if (!fArrTime)
    118117        return kFALSE;
    119 */
     118
    120119    return kTRUE;
    121120}
     
    123122
    124123// --------------------------------------------------------------------------
    125 // Evaluation of the mean arrival times (for now it stands for the maximum in slices units)
     124// Evaluation of the mean arrival times (spline interpolation)
    126125// per pixel and store them in the MArrivalTime container.
    127126//
     
    130129    MRawEvtPixelIter pixel(fRawEvt);
    131130
    132   //fArrTime->Calc((const MRawEvtData&) *fRawEvt,(const MGeomCam&) *fGeom);
    133   //fArrTime->SetReadyToSave();
     131    fArrTime->Calc((const MRawEvtData&) *fRawEvt,(const MGeomCam&) *fGeom);
     132    fArrTime->SetReadyToSave();
    134133
    135   //while (pixel.Next()){;}
     134    while (pixel.Next()){;}
    136135   
    137136    return kTRUE;
  • trunk/MagicSoft/Mars/manalysis/MArrivalTimeCalc.h

    r2659 r2752  
    66// MArrivalTimeCalc                                                        //
    77//                                                                         //
    8 // Evaluates the number of time slice into which the signal reaches a max. //
    9 // P R E L I M I N A R Y                                                   //
    10 // Other more sophisticated methods have to be implemented.                //
     8// Evaluates the Arrival Times                                            //
     9//                                                                         //
     10//                                                                         //
    1111/////////////////////////////////////////////////////////////////////////////
    1212
     
    2727    MRawRunHeader  *fRunHeader;  // RunHeader information
    2828    MGeomCam       *fGeom;       // Geometry information
    29 //  MArrivalTime   *fArrTime;    // Container with the photons arrival times
     29    MArrivalTime   *fArrTime;    // Container with the photons arrival times
    3030
    3131    Bool_t          fEnableFix;  // fix for a bug in files from older camera versions (<=40)
     
    4242    MArrivalTimeCalc(const char *name=NULL, const char *title=NULL);
    4343
    44     // FIXME: The array size should be checked!
    45 
    4644    ~MArrivalTimeCalc(){}
    4745
    48     ClassDef(MArrivalTimeCalc, 0)   // Task to calculate cerenkov photons from raw data
     46    ClassDef(MArrivalTimeCalc, 0)   // Task to calculate Arrival Times from raw data
    4947};
    5048
Note: See TracChangeset for help on using the changeset viewer.