EMCAL/AliCaloRawAnalyzer.cxx

   1 /**************************************************************************
   2  * This file is property of and copyright by                              *
   3  * the Relativistic Heavy Ion Group (RHIG), Yale University, US, 2009     *
   4  *                                                                        *
   5  * Primary Author: Per Thomas Hille <perthomas.hille@yale.edu>            *
   6  *                                                                        *
   7  * Contributors are mentioned in the code where appropriate.              *
   8  * Please report bugs to p.t.hille@fys.uio.no                             *
   9  *                                                                        *
  10  * Permission to use, copy, modify and distribute this software and its   *
  11  * documentation strictly for non-commercial purposes is hereby granted   *
  12  * without fee, provided that the above copyright notice appears in all   *
  13  * copies and that both the copyright notice and this permission notice   *
  14  * appear in the supporting documentation. The authors make no claims     *
  15  * about the suitability of this software for any purpose. It is          *
  16  * provided "as is" without express or implied warranty.                  *
  17  **************************************************************************/
  18
  19
  20 // Base class for extraction
  21 // of signal amplitude and peak position
  22 // From CALO Calorimeter RAW data (from the RCU)
  23 // Contains some utilities for preparing / selecting
  24 // Signals suitable for signal extraction
  25 // By derived classes
  26
  27 #include "AliLog.h"
  28 #include "AliCaloRawAnalyzer.h"
  29 #include "AliCaloBunchInfo.h"
  30 #include "AliCaloFitResults.h"
  31 #include "TMath.h"
  32 #include <iostream>
  33 using namespace std;
  34
  35 //ClassImp(AliCaloRawAnalyzer)
  36
  37 AliCaloRawAnalyzer::AliCaloRawAnalyzer(const char *name, const char *nameshort) :  TObject(),
  38                                                                                    fMinTimeIndex(-1),
  39                                                                                    fMaxTimeIndex(-1),
  40                                                                                    fFitArrayCut(5),
  41                                                                                    fAmpCut(4),
  42                                                                                    fNsampleCut(5),
  43                                                                                    fIsZerosupressed( false ),
  44                                                                                    fVerbose( false )
  45 {
  46   //Comment
  47   sprintf(fName, "%s", name);
  48   sprintf(fNameShort, "%s", nameshort);
  49
  50   for(int i=0; i < MAXSAMPLES; i++ )
  51     {
  52       fReversed[i] = 0;
  53     }
  54 }
  55
  56 AliCaloRawAnalyzer::~AliCaloRawAnalyzer()
  57 {
  58
  59 }
  60
  61
  62 void
  63 AliCaloRawAnalyzer::SetTimeConstraint(const int min, const int max )
  64 {
  65   //Require that the bin if the maximum ADC value is between min and max (timebin)
  66   if(  ( min > max ) || min > MAXSAMPLES  || max > MAXSAMPLES  )
  67     {
  68       AliWarning( Form( "Attempt to set Invalid time bin range (Min , Max) = (%d, %d), Ingored",  min, max ) );
  69     }
  70   else
  71     {
  72       fMinTimeIndex  = min;
  73       fMaxTimeIndex  = max;
  74     }
  75 }
  76
  77
  78 UShort_t
  79 AliCaloRawAnalyzer::Max(const UShort_t *data, const int length ) const
  80 {
  81   //------------
  82   UShort_t tmpmax  = data[0];
  83
  84   for(int i=0; i < length; i++)
  85     {
  86       if( tmpmax  <  data[i] )
  87         {
  88           tmpmax = data[i];
  89         }
  90     }
  91   return tmpmax;
  92 }
  93
  94
  95 void
  96 AliCaloRawAnalyzer::SelectSubarray( const Double_t *fData, const int length, const short maxindex, int *const first,  int *const last ) const
  97 {
  98   //Selection of subset of data from one bunch that will be used for fitting or
  99   //Peak finding.  Go to the left and right of index of the maximum time bin
 100   //Untile the ADC value is less that fFitArrayCut
 101   int tmpfirst =  maxindex;
 102   int tmplast  =  maxindex;
 103
 104   while((  tmpfirst  ) > 0  &&  ( fData[tmpfirst] >  fFitArrayCut   ))
 105     {
 106       tmpfirst -- ;
 107     }
 108
 109   while(( tmplast ) <  length   && ( fData [tmplast] >  fFitArrayCut ))
 110     {
 111       tmplast ++;
 112     }
 113
 114   *first = tmpfirst +1;
 115   *last =  tmplast -1;
 116 }
 117
 118
 119
 120 Float_t
 121 AliCaloRawAnalyzer::ReverseAndSubtractPed( const AliCaloBunchInfo *bunch, const UInt_t /*altrocfg1*/,  const UInt_t /*altrocfg2*/, double *outarray ) const
 122 {
 123   //Time sample comes in reversed order, revers them back
 124   //Subtract the baseline based on content of altrocfg1 and altrocfg2.
 125   Int_t length =  bunch->GetLength();
 126   const UShort_t *sig = bunch->GetData();
 127
 128   double ped = EvaluatePedestal( sig , length);
 129
 130   for( int i=0; i < length; i++ )
 131     {
 132       outarray[i] = sig[length -i -1] - ped;
 133     }
 134
 135   return ped;
 136 }
 137
 138
 139
 140 Float_t
 141 AliCaloRawAnalyzer::EvaluatePedestal(const UShort_t * const data, const int /*length*/ ) const
 142 {
 143   //  double ped = 0;
 144   double tmp = 0;
 145
 146   if( fIsZerosupressed == false )
 147     {
 148       for(int i=0; i < 5; i++ )
 149         {
 150           tmp += data[i];
 151         }
 152     }
 153
 154   //  cout << __FILE__ << __LINE__ << "XXXXXXXXXXX returning " <<   tmp/5 << endl;
 155
 156   return  tmp/5;
 157 }
 158
 159
 160 short
 161 AliCaloRawAnalyzer::Max( const AliCaloBunchInfo *const bunch , int *const maxindex ) const
 162 {
 163   //comment
 164   short tmpmax = -1;
 165   int tmpindex = -1;
 166   const UShort_t *sig = bunch->GetData();
 167
 168   for(int i=0; i < bunch->GetLength(); i++ )
 169     {
 170       if( sig[i] > tmpmax )
 171         {
 172           tmpmax   =  sig[i];
 173           tmpindex =  i;
 174         }
 175     }
 176
 177   if(maxindex != 0 )
 178     {
 179       //   *maxindex =  bunch->GetLength() -1 - tmpindex + bunch->GetStartBin();
 180        *maxindex =  bunch->GetLength() -1 - tmpindex + bunch->GetStartBin();
 181     }
 182
 183   return  tmpmax;
 184 }
 185
 186
 187 int
 188 AliCaloRawAnalyzer::SelectBunch( const vector<AliCaloBunchInfo> &bunchvector,short *const maxampbin, short *const maxamplitude ) const
 189 {
 190   //We select the bunch with the highest amplitude unless any time constraints is set
 191   short max = -1;
 192   short bunchindex = -1;
 193   short maxall = -1;
 194   int indx = -1;
 195
 196   for(unsigned int i=0; i < bunchvector.size(); i++ )
 197     {
 198       max = Max(  &bunchvector.at(i), &indx ); // CRAP PTH, bug fix, trouble if more than one bunches
 199       if( IsInTimeRange( indx) )
 200         {
 201           if( max > maxall )
 202             {
 203               maxall = max;
 204               bunchindex = i;
 205               *maxampbin     = indx;
 206               *maxamplitude  = max;
 207             }
 208         }
 209     }
 210
 211
 212   //  *maxampbin     = indx;
 213   //  *maxamplitude  = max;
 214
 215   return  bunchindex;
 216 }
 217
 218
 219 bool
 220 AliCaloRawAnalyzer::IsInTimeRange( const int maxindex  ) const
 221 {
 222   // Ckeck if the index of the max ADC vaue is consistent with trigger.
 223   if( ( fMinTimeIndex  < 0 && fMaxTimeIndex  < 0) ||fMaxTimeIndex  < 0 )
 224     {
 225       return true;
 226     }
 227   return ( maxindex < fMaxTimeIndex ) && ( maxindex > fMinTimeIndex  ) ? true : false;
 228 }
 229
 230
 231 void
 232 AliCaloRawAnalyzer::PrintBunches( const vector<AliCaloBunchInfo> &bvctr ) const
 233 {
 234   //comment
 235   cout << __FILE__ << __LINE__<< "*************** Printing Bunches *******************" << endl;
 236   cout << __FILE__ << __LINE__<< "*** There are " << bvctr.size() << ", bunches" << endl;
 237
 238   for(unsigned int i=0; i < bvctr.size() ; i++ )
 239     {
 240       PrintBunch( bvctr.at(i) );
 241       cout << " bunch = "  <<  i  << endl;
 242     }
 243   cout << __FILE__ << __LINE__<< "*************** Done ... *******************" << endl;
 244 }
 245
 246
 247 void
 248 AliCaloRawAnalyzer::PrintBunch( const AliCaloBunchInfo &bunch ) const
 249 {
 250   //comment
 251   cout << __FILE__ << ":" << __LINE__ << endl;
 252   cout << __FILE__ << __LINE__   << ", startimebin = " << bunch.GetStartBin() << ", length =" <<  bunch.GetLength() << endl;
 253   const UShort_t *sig =  bunch.GetData();
 254
 255   for ( Int_t j = 0; j <  bunch.GetLength();  j++)
 256     {
 257       printf("%d\t", sig[j] );
 258     }
 259   cout << endl;
 260 }
 261
 262
 263 AliCaloFitResults
 264 AliCaloRawAnalyzer::Evaluate( const vector<AliCaloBunchInfo>  &/*bunchvector*/, const UInt_t /*altrocfg1*/,  const UInt_t /*altrocfg2*/)
 265 { // method to do the selection of what should possibly be fitted
 266   // not implemented for base class
 267   return AliCaloFitResults( 0, 0, 0, 0, 0, 0, 0 );
 268 }
 269
 270
 271 int
 272 AliCaloRawAnalyzer::PreFitEvaluateSamples( const vector<AliCaloBunchInfo>  &bunchvector, const UInt_t altrocfg1,  const UInt_t altrocfg2, Int_t & index, Float_t & maxf, short & maxamp, short & maxampindex, Float_t & ped, int & first, int & last)
 273 { // method to do the selection of what should possibly be fitted
 274   int nsamples  = 0;
 275   index = SelectBunch( bunchvector,  &maxampindex,  &maxamp ); // select the bunch with the highest amplitude unless any time constraints is set
 276
 277
 278   if( index >= 0 && maxamp > fAmpCut) // something valid was found, and non-zero amplitude
 279     {
 280       // use more convenient numbering and possibly subtract pedestal
 281       ped  = ReverseAndSubtractPed( &(bunchvector.at(index)),  altrocfg1, altrocfg2, fReversed  );
 282       maxf = TMath::MaxElement( bunchvector.at(index).GetLength(),  fReversed );
 283
 284       if ( maxf > fAmpCut ) // possibly significant signal
 285         {
 286           // select array around max to possibly be used in fit
 287           maxampindex -= bunchvector.at(index).GetStartBin(); // PTH - why isn't this index also reversed for call below?
 288           SelectSubarray( fReversed,  bunchvector.at(index).GetLength(),  maxampindex , &first, &last);
 289
 290           // sanity check: maximum should not be in first or last bin
 291           // if we should do a fit
 292           if (first!=maxampindex && last!=maxampindex) {
 293             // calculate how many samples we have
 294             nsamples =  last - first + 1;
 295           }
 296         }
 297     }
 298
 299   return nsamples;
 300 }
 301