Ignore:
Timestamp:
08/16/04 15:26:49 (20 years ago)
Author:
tbretz
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars/mpedestal
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.cc

    r4613 r4629  
    394394
    395395  const Int_t npixels  = fPedestals->GetSize();
     396  const Int_t areas    = fPedestals->GetAverageAreas();
     397  const Int_t sectors  = fPedestals->GetAverageSectors();
    396398 
    397399  if (fSumx.GetSize()==0)
    398   {
     400    {
    399401      fSumx. Set(npixels);
    400402      fSumx2.Set(npixels);
    401 
     403     
     404      fAreaSumx. Set(areas);
     405      fAreaSumx2.Set(areas);
     406      fAreaValid.Set(areas);
     407     
     408      fSectorSumx. Set(sectors);
     409      fSectorSumx2.Set(sectors);
     410      fSectorValid.Set(sectors);
     411     
    402412      fSumx.Reset();
    403413      fSumx2.Reset();
    404   }
     414    }
    405415 
    406416  if (fWindowSizeHiGain == 0 && fWindowSizeLoGain == 0)
    407417    {
    408         *fLog << err << GetDescriptor()
     418      *fLog << err << GetDescriptor()
    409419            << ": Number of extracted Slices is 0, abort ... " << endl;
    410         return kFALSE;
     420      return kFALSE;
    411421    }
    412422 
     
    429439Int_t MPedCalcPedRun::Process()
    430440{
    431     MRawEvtPixelIter pixel(fRawEvt);
    432 
    433     while (pixel.Next())
    434     {
    435         const UInt_t idx    = pixel.GetPixelId();
    436 
    437         Byte_t *ptr = pixel.GetHiGainSamples() + fHiGainFirst;
    438         Byte_t *end = ptr + fWindowSizeHiGain;
    439 
    440         UInt_t sum = 0;
    441         UInt_t sqr = 0;
    442 
    443         if (fWindowSizeHiGain != 0)
    444         {
    445             do
    446             {
    447                 sum += *ptr;
    448                 sqr += *ptr * *ptr;
    449             }
    450             while (++ptr != end);
    451         }
    452 
    453         if (fWindowSizeLoGain != 0)
    454         {
    455             ptr = pixel.GetLoGainSamples() + fLoGainFirst;
    456             end = ptr + fWindowSizeLoGain;
    457 
    458             do
    459             {
    460                 sum += *ptr;
    461                 sqr += *ptr * *ptr;
    462             }
    463             while (++ptr != end);
    464         }
    465 
    466         const Float_t msum = (Float_t)sum;
    467 
    468         //
    469         // These three lines have been uncommented by Markus Gaug
    470         // If anybody needs them, please contact me!!
    471         //
    472         //      const Float_t higainped = msum/fNumHiGainSlices;
    473         //      const Float_t higainrms = TMath::Sqrt((msqr-msum*msum/fNumHiGainSlices)/(fNumHiGainSlices-1.));
    474         //      (*fPedestals)[idx].Set(higainped, higainrms);
    475 
    476         fSumx[idx]          += msum;
    477         //
    478         // The old version:
    479         //
    480         //       const Float_t msqr = (Float_t)sqr;
    481         //      fSumx2[idx] += msqr;
    482         //
    483         // The new version:
    484         //
    485         const Float_t sqrsum  = msum*msum;
    486         fSumx2[idx]          += sqrsum;
    487     }
    488 
    489     fNumSamplesTot += fWindowSizeHiGain + fWindowSizeLoGain;
    490 
    491     return kTRUE;
    492 /*
    493441  MRawEvtPixelIter pixel(fRawEvt);
    494442 
     
    497445      const UInt_t idx    = pixel.GetPixelId();
    498446      const UInt_t aidx   = (*fGeom)[idx].GetAidx();
    499       const UInt_t sector = (*fGeom)[idx].GetSector();
     447      const UInt_t sector = (*fGeom)[idx].GetSector();     
    500448
    501449      Byte_t *ptr = pixel.GetHiGainSamples() + fHiGainFirst;
     
    542490      fSumx[idx]          += msum;
    543491      fAreaSumx[aidx]     += msum;
    544       fSectorSumx[sector] += msum;
     492      fSectorSumx[sector] += msum;     
    545493      //
    546494      // The old version:
     
    554502      fSumx2[idx]          += sqrsum;
    555503      fAreaSumx2[aidx]     += sqrsum;
    556       fSectorSumx2[sector] += sqrsum;
     504      fSectorSumx2[sector] += sqrsum;     
    557505    }
    558506 
    559507  fPedestals->SetReadyToSave();
    560508  fNumSamplesTot += fWindowSizeHiGain + fWindowSizeLoGain;
    561 */
     509
    562510  return kTRUE;
    563511}
     
    569517Int_t MPedCalcPedRun::PostProcess()
    570518{
    571     const ULong_t nevts = GetNumExecutions();
    572     if (nevts<1)
    573     {
    574         *fLog << err << "ERROR - Not enough events recorded...abort." << endl;
    575         return kFALSE;
    576     }
    577 
    578     *fLog << flush << inf << "Calculating pedestals..." << flush;
    579 
    580     // Compute pedestals and rms from the whole run
    581     const ULong_t n     = fNumSamplesTot;
    582     const ULong_t npix  = fGeom->GetNumPixels();
    583 
    584     for (UInt_t pixidx=0; pixidx<npix; pixidx++)
    585     {
    586         const Float_t sum  = fSumx.At(pixidx);
    587         const Float_t sum2 = fSumx2.At(pixidx);
    588         const Float_t higainped = sum/n;
    589         //
    590         // The old version:
    591         //
    592         //      const Float_t higainrms = TMath::Sqrt((sum2-sum*sum/n)/(n-1.));
    593         //
    594         // The new version:
    595         //
    596         // 1. Calculate the Variance of the sums:
    597         Float_t higainVar = (sum2-sum*sum/nevts)/(nevts-1.);
    598         // 2. Scale the variance to the number of slices:
    599         higainVar /= (Float_t)(fWindowSizeHiGain+fWindowSizeLoGain);
    600         // 3. Calculate the RMS from the Variance:
    601         (*fPedestals)[pixidx].Set(higainped, higainVar < 0 ? 0. : TMath::Sqrt(higainVar), 0, nevts);
    602     }
    603 
    604     *fLog << flush << inf << "Calculating means..." << flush;
    605 
    606     fPedestals->SetTotalEntries(fNumSamplesTot);
    607     fPedestals->ReCalc(*fGeom);
    608     fPedestals->SetReadyToSave();
    609 
    610     return kTRUE;
     519  // Compute pedestals and rms from the whole run
     520  const ULong_t n     = fNumSamplesTot;
     521  const ULong_t nevts = GetNumExecutions();
     522
     523  MRawEvtPixelIter pixel(fRawEvt);
     524 
     525  while (pixel.Next())
     526    {
     527
     528      const Int_t  pixid  = pixel.GetPixelId();
     529      const UInt_t aidx   = (*fGeom)[pixid].GetAidx();
     530      const UInt_t sector = (*fGeom)[pixid].GetSector();     
     531     
     532      fAreaValid  [aidx]++;
     533      fSectorValid[sector]++;
     534
     535      const Float_t sum  = fSumx.At(pixid);
     536      const Float_t sum2 = fSumx2.At(pixid);
     537      const Float_t higainped = sum/n;
     538      //
     539      // The old version:
     540      //
     541      //      const Float_t higainrms = TMath::Sqrt((sum2-sum*sum/n)/(n-1.));
     542      //
     543      // The new version:
     544      //
     545      // 1. Calculate the Variance of the sums:
     546      Float_t higainVar = (sum2-sum*sum/nevts)/(nevts-1.);
     547      // 2. Scale the variance to the number of slices:
     548      higainVar /= (Float_t)(fWindowSizeHiGain+fWindowSizeLoGain);
     549      // 3. Calculate the RMS from the Variance:
     550      (*fPedestals)[pixid].Set(higainped, higainVar < 0 ? 0. : TMath::Sqrt(higainVar));
     551
     552    }
    611553
    612554  //
    613555  // Loop over the (two) area indices to get the averaged pedestal per aidx
    614556  //
    615   /*
    616 
    617    // Compute pedestals and rms from the whole run
    618     const ULong_t n     = fNumSamplesTot;
    619     const ULong_t nevts = GetNumExecutions();
    620 
    621     MRawEvtPixelIter pixel(fRawEvt);
    622 
    623     while (pixel.Next())
    624     {
    625         const Int_t  pixid  = pixel.GetPixelId();
    626         const UInt_t aidx   = (*fGeom)[pixid].GetAidx();
    627         const UInt_t sector = (*fGeom)[pixid].GetSector();
    628 
    629         fAreaValid  [aidx]++;
    630         fSectorValid[sector]++;
    631 
    632         const Float_t sum  = fSumx.At(pixid);
    633         const Float_t sum2 = fSumx2.At(pixid);
    634         const Float_t higainped = sum/n;
    635         //
    636         // The old version:
    637         //
    638         //      const Float_t higainrms = TMath::Sqrt((sum2-sum*sum/n)/(n-1.));
    639         //
    640         // The new version:
    641         //
    642         // 1. Calculate the Variance of the sums:
    643         Float_t higainVar = (sum2-sum*sum/nevts)/(nevts-1.);
    644         // 2. Scale the variance to the number of slices:
    645         higainVar /= (Float_t)(fWindowSizeHiGain+fWindowSizeLoGain);
    646         // 3. Calculate the RMS from the Variance:
    647         (*fPedestals)[pixid].Set(higainped, higainVar < 0 ? 0. : TMath::Sqrt(higainVar), 0, nevts);
    648     }
    649557  for (Int_t aidx=0; aidx<fAreaValid.GetSize(); aidx++)
    650558    {
     
    674582      fPedestals->GetAverageArea(aidx).Set(higainped, higainrms);
    675583    }
    676     */
     584 
    677585  //
    678586  // Loop over the (six) sector indices to get the averaged pedestal per sector
    679587  //
    680   /*
    681588  for (Int_t sector=0; sector<fSectorValid.GetSize(); sector++)
    682589    {
     
    706613      fPedestals->GetAverageSector(sector).Set(higainped, higainrms);
    707614    }
     615 
    708616  fPedestals->SetTotalEntries(fNumSamplesTot);
    709617  fPedestals->SetReadyToSave();
    710618
    711619  return kTRUE;
    712     */
    713620}
    714621
  • trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.h

    r4609 r4629  
    3333    TArrayD fSumx;             // sum of values
    3434    TArrayD fSumx2;            // sum of squared values
    35     /*
    3635    TArrayD fAreaSumx;         // averaged sum of values per area idx
    3736    TArrayD fAreaSumx2;        // averaged sum of squared values per area idx
     
    4039    TArrayD fSectorSumx2;      // averaged sum of squared values per sector
    4140    TArrayI fSectorValid;      // number of valid pixel with sector idx
    42     */
    4341
    4442    Int_t  PreProcess (MParList *pList);
Note: See TracChangeset for help on using the changeset viewer.