Changeset 5558 for trunk/MagicSoft/Mars


Ignore:
Timestamp:
12/03/04 20:17:06 (20 years ago)
Author:
tbretz
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars
Files:
13 edited

Legend:

Unmodified
Added
Removed
  • trunk/MagicSoft/Mars/Changelog

    r5552 r5558  
    3030     - added possibility to set fPedestals directly by pointer
    3131     - unfolded some Form statements
     32
     33   * callisto.cc:
     34     - changed callisto to support second pedestal loop -- WORK IN PROGRESS
     35
     36   * mbase/MParContainer.[h,cc]:
     37     - moved code from GetDescriptor to a static function GetDescriptor
     38     
     39   * mbase/MStatusDisplay.cc:
     40     - removed a oboslete debug output
     41
     42   * mbase/MTaskEnv.h:
     43     - made ReadEnv public
     44
     45   * mcalib/MCalibrationChargeCam.[h,cc]:
     46     - changed some returned TArrayF* to TArrayF
     47
     48   * mimage/MImgCleanStd.cc:
     49     - changed point of output of name of MPedPhotCam
     50
     51   * mjobs/MJCalibration.cc:
     52     - don't use MParList to hand ExtractorCam over
     53     - some changes to the structur for Writing (only consmetics)
     54
     55   * mjobs/MJPedestal.[h,cc]:
     56     - ordered includes correctly
     57     - for the moment removed fNameExtractorCam -- simplifies thing a lot
     58     - make a clone of the extractor given by the user - for sanity
     59     - correct handling of the allocated memory for fExtractor
     60     - replaced counts/slice by cts/slice for space reasons
     61     - changed name of new tabs for space reasons
     62     - added possibility to CheckEnv to set extractor from resource-file
     63     - outsourced some code to SetupExtractor
     64     - given a name to MFillH
     65     - changed handling of SetNoiseCalculation which is now set correctly all
     66       the time
     67     - Writing must still be checked!!!
     68
     69   * mjobs/MJob.[h,cc]:
     70     - allow WriteContainer to use any TObject
     71     - added Getter-functions for the TEnv stuff
     72
     73   * mpedestal/MExtractPedestal.[h,cc]:
     74     - allow setting of fPedestalIn by pointer directly
     75     - fixed handling of fPedestalIn accordingly
     76     - changed some arguments from pointer to reference
     77     - changed some accesses to TArrays from At to []-operator
     78     - shortened and enhanced output
     79     - fixed place and type of screen output
     80
     81   * mpedestal/MPedCalcFromLoGain.[h,cc]:
     82     - removed GetSlice -> replaced by a array with a copy of the data
     83     - changed some loops to pointer arithmetic for speed reasons
     84       in this case
     85
     86   * mpedestal/MPedCalcPedRun.[h,cc]:
     87     - changed handling of first pedestal run for simplicity
     88     - replaced some Form calls
     89     - changed direct handling of single bits to correct enums
     90     - fixed wrong 'all' in output
     91     - simplified output
     92     - fixed some stuff in the output
     93
     94   * mpedestal/MPedestalCam.[h,cc]:
     95     - fixed the Copy function - WITHOUT this fix the calibration could
     96       not have worked properly at all.
     97     - changed some function names to Mars standrad names
     98     - removed some obsolete loops - replaced by ForEach
     99     - changed some TArrayF* returnes to TArrayF
     100
     101   * msignal/MExtractTime.cc, msignal/MExtractTimeAndCharge.cc,
     102     msignal/MExtractTimeAndChargeDigitalFilter.cc,
     103     msignal/MExtractor.cc
     104     - fixed some wrong debug levels in output
     105     - some simplification and shortening to output
     106
    32107
    33108
  • trunk/MagicSoft/Mars/mpedestal/MExtractPedestal.cc

    r5521 r5558  
    150150/////////////////////////////////////////////////////////////////////////////
    151151#include "MExtractPedestal.h"
    152 #include "MExtractTimeAndCharge.h"
    153152
    154153#include "MParList.h"
     
    168167#include "MGeomCam.h"
    169168
     169#include "MTaskEnv.h"
     170#include "MExtractTimeAndCharge.h"
     171
    170172ClassImp(MExtractPedestal);
    171173
     
    173175
    174176const TString MExtractPedestal::fgNamePedestalCam = "MPedestalCam";
     177
    175178// --------------------------------------------------------------------------
    176179//
     
    189192      fExtractWinFirst(0), fExtractWinSize(0)
    190193{
    191 
    192   fName  = name  ? name  : "MExtractPedestal";
    193   fTitle = title ? title : "Base class to calculate pedestals";
    194  
    195   AddToBranchList("fHiGainPixId");
    196   AddToBranchList("fLoGainPixId");
    197   AddToBranchList("fHiGainFadcSamples");
    198   AddToBranchList("fLoGainFadcSamples");
    199  
    200   SetPedestalUpdate(kTRUE);
    201 
    202   SetNamePedestalCamIn();
    203   SetNamePedestalCamOut();
    204   SetNumEventsDump();
    205   SetNumAreasDump();
    206   SetNumSectorsDump();
    207  
    208   Clear();
     194    fName  = name  ? name  : "MExtractPedestal";
     195    fTitle = title ? title : "Base class to calculate pedestals";
     196
     197    AddToBranchList("fHiGainPixId");
     198    AddToBranchList("fLoGainPixId");
     199    AddToBranchList("fHiGainFadcSamples");
     200    AddToBranchList("fLoGainFadcSamples");
     201
     202    SetPedestalUpdate(kTRUE);
     203
     204    SetNamePedestalCamIn();
     205    SetNamePedestalCamOut();
     206    SetNumEventsDump();
     207    SetNumAreasDump();
     208    SetNumSectorsDump();
     209
     210    Clear();
    209211}
    210212
     
    212214{
    213215    // Reset contents of arrays.
    214   fSumx.Reset();       
    215   fSumx2.Reset();       
    216   fSumAB0.Reset();     
    217   fSumAB1.Reset();     
    218   fAreaSumx.Reset();   
    219   fAreaSumx2.Reset();   
    220   fAreaSumAB0.Reset(); 
    221   fAreaSumAB1.Reset(); 
    222   fAreaFilled.Reset();   
    223   fAreaValid.Reset();   
    224   fSectorSumx.Reset(); 
    225   fSectorSumx2.Reset();
    226   fSectorSumAB0.Reset();
    227   fSectorSumAB1.Reset();
    228   fSectorFilled.Reset();
    229   fSectorValid.Reset();
     216    fSumx.Reset();
     217    fSumx2.Reset();
     218    fSumAB0.Reset();
     219    fSumAB1.Reset();
     220    fAreaSumx.Reset();
     221    fAreaSumx2.Reset();
     222    fAreaSumAB0.Reset();
     223    fAreaSumAB1.Reset();
     224    fAreaFilled.Reset();
     225    fAreaValid.Reset();
     226    fSectorSumx.Reset();
     227    fSectorSumx2.Reset();
     228    fSectorSumAB0.Reset();
     229    fSectorSumAB1.Reset();
     230    fSectorFilled.Reset();
     231    fSectorValid.Reset();
    230232
    231233}
     
    248250  fRunHeader    = NULL;
    249251  fEvtHeader    = NULL;
    250   fPedestalsIn  = NULL;
    251252  fPedestalsOut = NULL;
    252253
     
    270271 
    271272  if (odd)
    272     {
     273  {
    273274      *fLog << warn << GetDescriptor();
    274       *fLog << " - WARNING: Window size in SetExtraxtWindow has to be even... " 
    275         " raising from " << windows << " to " << flush;
     275      *fLog << " - WARNING: Window size in SetExtraxtWindow has to be even... ";
     276      *fLog << " raising from " << windows << " to ";
    276277      windows += 1;
    277278      *fLog << windows << "!" << endl;
    278279      rc = kFALSE;
    279     }
    280 
    281   if (windows==0) 
    282     {
     280  }
     281
     282  if (windows==0)
     283  {
    283284      *fLog << warn << GetDescriptor();
    284285      *fLog << " - WARNING: Window size in SetExtraxtWindow has to be > 0... adjusting to 2!" << endl;
    285286      windows = 2;
    286287      rc = kFALSE;
    287     }
     288  }
    288289
    289290  fExtractWinSize  = windows;
     
    315316  fRawEvt = (MRawEvtData*)pList->FindObject(AddSerialNumber("MRawEvtData"));
    316317  if (!fRawEvt)
    317     {
     318  {
    318319      *fLog << err << AddSerialNumber("MRawEvtData") << " not found... aborting." << endl;
    319320      return kFALSE;
    320     }
     321  }
    321322 
    322323  fRunHeader = (MRawRunHeader*)pList->FindObject(AddSerialNumber("MRawRunHeader"));
    323324  if (!fRunHeader)
    324     {
     325  {
    325326      *fLog << err << AddSerialNumber("MRawRunHeader") << " not found... aborting." << endl;
    326327      return kFALSE;
    327     }
     328  }
    328329 
    329330  fEvtHeader = (MRawEvtHeader*)pList->FindObject(AddSerialNumber("MRawEvtHeader"));
    330331  if (!fEvtHeader)
    331     {
     332  {
    332333      *fLog << err << AddSerialNumber("MRawEvtHeader") << " not found... aborting." << endl;
    333334      return kFALSE;
    334     }
     335  }
    335336 
    336337  fGeom = (MGeomCam*)pList->FindObject(AddSerialNumber("MGeomCam"));
    337338  if (!fGeom)
    338     {
     339  {
    339340      *fLog << err << AddSerialNumber("MGeomCam") << " not found... aborting." << endl;
    340341      return kFALSE;
    341     }
    342  
    343   fPedestalsIn = (MPedestalCam*)pList->FindCreateObj("MPedestalCam", AddSerialNumber(fNamePedestalCamIn.Data()));
    344   if (!fPedestalsIn)
    345     {
    346       *fLog << err << fNamePedestalCamIn.Data() << " could not be found nor created... aborting" << endl;
     342  }
     343
     344  if (fExtractor && !fPedestalsIn)
     345  {
     346      fPedestalsIn = (MPedestalCam*)pList->FindObject("MPedestalCam", AddSerialNumber(fNamePedestalCamIn));
     347      if (!fPedestalsIn)
     348      {
     349          *fLog << err << AddSerialNumber(fNamePedestalCamIn) << " not found... aborting." << endl;
     350          return kFALSE;
     351      }
     352  }
     353
     354  fPedestalsOut = (MPedestalCam*)pList->FindCreateObj("MPedestalCam", AddSerialNumber(fNamePedestalCamOut));
     355  if (!fPedestalsOut)
    347356      return kFALSE;
    348     }
    349  
    350   fPedestalsOut = (MPedestalCam*)pList->FindCreateObj("MPedestalCam", AddSerialNumber(fNamePedestalCamOut.Data()));
    351   if (!fPedestalsOut)
    352     {
    353       *fLog << err << fNamePedestalCamOut.Data() << " could not be found nor created... aborting" << endl;
    354       return kFALSE;
    355     }
    356  
     357
    357358  *fLog << inf;
    358  
     359  Print();
     360
    359361  return kTRUE;
    360362}
     
    382384Bool_t MExtractPedestal::ReInit(MParList *pList)
    383385{
    384 
    385386  // If the size is not yet set, set the size
    386387  if (fSumx.GetSize()==0)
    387388  {
    388     const Int_t npixels  = fPedestalsOut->GetSize();
    389     const Int_t areas    = fPedestalsOut->GetAverageAreas();
    390     const Int_t sectors  = fPedestalsOut->GetAverageSectors();
    391    
    392     fSumx.  Set(npixels);
    393     fSumx2. Set(npixels);
    394     fSumAB0.Set(npixels);
    395     fSumAB1.Set(npixels);
    396    
    397     fAreaSumx.  Set(areas);
    398     fAreaSumx2. Set(areas);
    399     fAreaSumAB0.Set(areas);
    400     fAreaSumAB1.Set(areas);
    401     fAreaFilled.Set(areas);
    402     fAreaValid .Set(areas);
    403    
    404     fSectorSumx.  Set(sectors);
    405     fSectorSumx2. Set(sectors);
    406     fSectorSumAB0.Set(sectors);
    407     fSectorSumAB1.Set(sectors);
    408     fSectorFilled.Set(sectors);
    409     fSectorValid .Set(sectors);
    410 
    411     for (Int_t i=0; i<npixels; i++)
     389      const Int_t npixels  = fPedestalsOut->GetSize();
     390      const Int_t areas    = fPedestalsOut->GetNumAverageArea();
     391      const Int_t sectors  = fPedestalsOut->GetNumAverageSector();
     392
     393      fSumx.  Set(npixels);
     394      fSumx2. Set(npixels);
     395      fSumAB0.Set(npixels);
     396      fSumAB1.Set(npixels);
     397
     398      fAreaSumx.  Set(areas);
     399      fAreaSumx2. Set(areas);
     400      fAreaSumAB0.Set(areas);
     401      fAreaSumAB1.Set(areas);
     402      fAreaFilled.Set(areas);
     403      fAreaValid .Set(areas);
     404
     405      fSectorSumx.  Set(sectors);
     406      fSectorSumx2. Set(sectors);
     407      fSectorSumAB0.Set(sectors);
     408      fSectorSumAB1.Set(sectors);
     409      fSectorFilled.Set(sectors);
     410      fSectorValid .Set(sectors);
     411
     412      for (Int_t i=0; i<npixels; i++)
    412413      {
    413         const UInt_t aidx   = (*fGeom)[i].GetAidx();
    414         const UInt_t sector = (*fGeom)[i].GetSector();     
    415        
    416         fAreaValid  [aidx]  ++;
    417         fSectorValid[sector]++;
     414          const UInt_t aidx   = (*fGeom)[i].GetAidx();
     415          const UInt_t sector = (*fGeom)[i].GetSector();
     416
     417          fAreaValid  [aidx]  ++;
     418          fSectorValid[sector]++;
    418419      }
    419420  }
    420  
     421
    421422  if (fExtractor)
    422     {
     423  {
    423424      if (!fExtractor->InitArrays())
    424         return kFALSE;
    425       SetExtractWindow(fExtractor->GetHiGainFirst(),(Int_t)fExtractor->GetNumHiGainSamples());
    426     }
    427 
    428   Print();
     425          return kFALSE;
     426
     427      SetExtractWindow(fExtractor->GetHiGainFirst(), (Int_t)TMath::Nint(fExtractor->GetNumHiGainSamples()));
     428  }
    429429
    430430  return kTRUE;
    431431}
     432
     433Int_t MExtractPedestal::PostProcess()
     434{
     435    fPedestalsIn = NULL;
     436    return kTRUE;
     437}
     438
    432439
    433440// --------------------------------------------------------------------------
     
    441448Int_t MExtractPedestal::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
    442449{
    443 
    444   Bool_t rc=kFALSE;
    445  
    446   // find resource for numeventsdump
    447   if (IsEnvDefined(env, prefix, "NumEventsDump", print))
    448     {
    449       SetNumEventsDump(GetEnvValue(env, prefix, "NumEventsDump", (Int_t)fNumEventsDump));
    450       rc = kTRUE;
    451     }
    452 
    453   // find resource for numeventsdump
    454   if (IsEnvDefined(env, prefix, "NumAreasDump", print))
    455     {
    456       SetNumAreasDump(GetEnvValue(env, prefix, "NumAreasDump", (Int_t)fNumAreasDump));
    457       rc = kTRUE;
    458     }
    459 
    460   // find resource for numeventsdump
    461   if (IsEnvDefined(env, prefix, "NumSectorsDump", print))
    462     {
    463       SetNumSectorsDump(GetEnvValue(env, prefix, "NumSectorsDump", (Int_t)fNumSectorsDump));
    464       rc = kTRUE;
    465     }
    466 
    467   // find resource for pedestal update
    468   if (IsEnvDefined(env, prefix, "PedestalUpdate", print))
    469     {
    470       SetPedestalUpdate(GetEnvValue(env, prefix, "PedestalUpdate", fPedestalUpdate));
    471       rc = kTRUE;
    472     }
    473  
    474   // Find resources for ExtractWindow
    475   Int_t ef = fExtractWinFirst;
    476   Int_t es = fExtractWinSize;
    477   if (IsEnvDefined(env, prefix, "ExtractWinFirst", print))
    478     {
    479       ef = GetEnvValue(env, prefix, "ExtractWinFirst", ef);
    480       rc = kTRUE;
    481     }
    482   if (IsEnvDefined(env, prefix, "ExtractWinSize", print))
    483     {
    484       es = GetEnvValue(env, prefix, "ExtractWinSize", es);
    485       rc = kTRUE;
    486     }
    487  
    488   SetExtractWindow(ef,es);
    489  
    490   // find resource for MPedestalCam
    491   if (IsEnvDefined(env, prefix, "NamePedestalCamIn", print))
    492     {
    493       SetNamePedestalCamIn(GetEnvValue(env, prefix, "NamePedestalCamIn", fNamePedestalCamIn));
    494       rc = kTRUE;
    495     }
    496  
    497   if (IsEnvDefined(env, prefix, "NamePedestalCamOut", print))
    498     {
    499       SetNamePedestalCamOut(GetEnvValue(env, prefix, "NamePedestalCamOut", fNamePedestalCamOut));
    500       rc = kTRUE;
    501     }
    502  
    503   return rc;
     450    Bool_t rc=kFALSE;
     451
     452    // find resource for numeventsdump
     453    if (IsEnvDefined(env, prefix, "NumEventsDump", print))
     454    {
     455        SetNumEventsDump(GetEnvValue(env, prefix, "NumEventsDump", (Int_t)fNumEventsDump));
     456        rc = kTRUE;
     457    }
     458
     459    // find resource for numeventsdump
     460    if (IsEnvDefined(env, prefix, "NumAreasDump", print))
     461    {
     462        SetNumAreasDump(GetEnvValue(env, prefix, "NumAreasDump", (Int_t)fNumAreasDump));
     463        rc = kTRUE;
     464    }
     465
     466    // find resource for numeventsdump
     467    if (IsEnvDefined(env, prefix, "NumSectorsDump", print))
     468    {
     469        SetNumSectorsDump(GetEnvValue(env, prefix, "NumSectorsDump", (Int_t)fNumSectorsDump));
     470        rc = kTRUE;
     471    }
     472
     473    // find resource for pedestal update
     474    if (IsEnvDefined(env, prefix, "PedestalUpdate", print))
     475    {
     476        SetPedestalUpdate(GetEnvValue(env, prefix, "PedestalUpdate", fPedestalUpdate));
     477        rc = kTRUE;
     478    }
     479
     480    // Find resources for ExtractWindow
     481    Int_t ef = fExtractWinFirst;
     482    Int_t es = fExtractWinSize;
     483    if (IsEnvDefined(env, prefix, "ExtractWinFirst", print))
     484    {
     485        ef = GetEnvValue(env, prefix, "ExtractWinFirst", ef);
     486        rc = kTRUE;
     487    }
     488    if (IsEnvDefined(env, prefix, "ExtractWinSize", print))
     489    {
     490        es = GetEnvValue(env, prefix, "ExtractWinSize", es);
     491        rc = kTRUE;
     492    }
     493
     494    SetExtractWindow(ef,es);
     495
     496    // find resource for MPedestalCam
     497    if (IsEnvDefined(env, prefix, "NamePedestalCamIn", print))
     498    {
     499        SetNamePedestalCamIn(GetEnvValue(env, prefix, "NamePedestalCamIn", fNamePedestalCamIn));
     500        rc = kTRUE;
     501    }
     502
     503    if (IsEnvDefined(env, prefix, "NamePedestalCamOut", print))
     504    {
     505        SetNamePedestalCamOut(GetEnvValue(env, prefix, "NamePedestalCamOut", fNamePedestalCamOut));
     506        rc = kTRUE;
     507    }
     508
     509    return rc;
    504510}
    505511
     
    516522void MExtractPedestal::CalcPixResults(const UInt_t nevts, const UInt_t pixid)
    517523{
    518  
    519   const Float_t sum         = fSumx.At(pixid);
    520   const Float_t sum2        = fSumx2.At(pixid);
    521 
    522   // 1. Calculate the mean of the sums:
    523   Float_t ped         = sum/nevts;
    524  
    525   // 2. Calculate the Variance of the sums:
    526   Float_t var = (sum2-sum*sum/nevts)/(nevts-1.);
    527 
    528   // 3. Calculate the amplitude of the 150MHz "AB" noise
    529   Float_t abOffs = (fSumAB0[pixid] - fSumAB1[pixid]) / nevts;
    530 
    531   // 4. Scale the mean, variance and AB-noise to the number of slices:
    532   ped    /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
    533   var    /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
    534   abOffs /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
    535  
    536   // 5. Calculate the RMS from the Variance:
    537   const Float_t rms = var<0 ? 0 : TMath::Sqrt(var);
    538  
    539   (*fPedestalsOut)[pixid].Set(ped, rms, abOffs, nevts);
     524    const Float_t sum  = fSumx[pixid];
     525    const Float_t sum2 = fSumx2[pixid];
     526
     527    // 1. Calculate the mean of the sums:
     528    Float_t ped        = sum/nevts;
     529
     530    // 2. Calculate the Variance of the sums:
     531    Float_t var = (sum2-sum*sum/nevts)/(nevts-1.);
     532
     533    // 3. Calculate the amplitude of the 150MHz "AB" noise
     534    Float_t abOffs = (fSumAB0[pixid] - fSumAB1[pixid]) / nevts;
     535
     536    // 4. Scale the mean, variance and AB-noise to the number of slices:
     537    ped    /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
     538    var    /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
     539    abOffs /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
     540
     541    // 5. Calculate the RMS from the Variance:
     542    const Float_t rms = var<0 ? 0 : TMath::Sqrt(var);
     543
     544    (*fPedestalsOut)[pixid].Set(ped, rms, abOffs, nevts);
    540545}
    541546
     
    552557void MExtractPedestal::CalcAreaResults(const UInt_t nevts, const UInt_t napix, const UInt_t aidx)
    553558{
    554  
    555   const Float_t sum         = fAreaSumx.At(aidx);
    556   const Float_t sum2        = fAreaSumx2.At(aidx);
    557 
    558   // 1. Calculate the mean of the sums:
    559   Float_t ped         = sum/nevts;     
    560   //
    561   // 2. Calculate the Variance of the sums:
    562   //
    563   Float_t var = (sum2-sum*sum/nevts)/(nevts-1.);
    564   //
    565   // 3. Calculate the amplitude of the 150MHz "AB" noise
    566   //
    567   Float_t abOffs = (fAreaSumAB0[aidx] - fAreaSumAB1[aidx]) / nevts;
    568   //
    569   // 4. Scale the mean, variance and AB-noise to the number of slices:
    570   //
    571   ped    /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
    572   var    /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
    573   abOffs /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
    574   //
    575   // 5. Scale the mean, variance and AB-noise to the number of pixels:
    576   //
    577   ped    /= napix;
    578   var    /= napix;
    579   abOffs /= napix;
    580   //
    581   // 6. Calculate the RMS from the Variance:
    582   //
    583   const Float_t rms = var<0 ? 0 : TMath::Sqrt(var);
    584 
    585   fPedestalsOut->GetAverageArea(aidx).Set(ped, rms,abOffs,nevts);
    586  
     559    const Float_t sum  = fAreaSumx[aidx];
     560    const Float_t sum2 = fAreaSumx2[aidx];
     561
     562    // 1. Calculate the mean of the sums:
     563    Float_t ped        = sum/nevts;
     564
     565    // 2. Calculate the Variance of the sums:
     566    Float_t var = (sum2-sum*sum/nevts)/(nevts-1.);
     567
     568    // 3. Calculate the amplitude of the 150MHz "AB" noise
     569    Float_t abOffs = (fAreaSumAB0[aidx] - fAreaSumAB1[aidx]) / nevts;
     570
     571    // 4. Scale the mean, variance and AB-noise to the number of slices:
     572    ped    /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
     573    var    /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
     574    abOffs /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
     575
     576    // 5. Scale the mean, variance and AB-noise to the number of pixels:
     577    ped    /= napix;
     578    var    /= napix;
     579    abOffs /= napix;
     580
     581    // 6. Calculate the RMS from the Variance:
     582    const Float_t rms = var<0 ? 0 : TMath::Sqrt(var);
     583
     584    fPedestalsOut->GetAverageArea(aidx).Set(ped, rms, abOffs, nevts);
    587585}
    588586
     
    599597void MExtractPedestal::CalcSectorResults(const UInt_t nevts, const UInt_t nspix, const UInt_t sector)
    600598{
    601  
    602   const Float_t sum         = fSectorSumx.At(sector);
    603   const Float_t sum2        = fSectorSumx2.At(sector);
    604 
    605   // 1. Calculate the mean of the sums:
    606   Float_t ped         = sum/nevts;     
    607   //
    608   // 2. Calculate the Variance of the sums:
    609   //
    610   Float_t var = (sum2-sum*sum/nevts)/(nevts-1.);
    611   //
    612   // 3. Calculate the amplitude of the 150MHz "AB" noise
    613   //
    614   Float_t abOffs = (fSectorSumAB0[sector] - fSectorSumAB1[sector]) / nevts;
    615   //
    616   // 4. Scale the mean, variance and AB-noise to the number of slices:
    617   //
    618   ped    /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
    619   var    /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
    620   abOffs /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
    621   //
    622   // 5. Scale the mean, variance and AB-noise to the number of pixels:
    623   //
    624   ped    /= nspix;
    625   var    /= nspix;
    626   abOffs /= nspix;
    627   //
    628   // 6. Calculate the RMS from the Variance:
    629   //
    630   const Float_t rms = var<0 ? 0 : TMath::Sqrt(var);
    631 
    632   fPedestalsOut->GetAverageSector(sector).Set(ped, rms,abOffs,nevts);
    633 }
    634 
    635 
     599    const Float_t sum  = fSectorSumx[sector];
     600    const Float_t sum2 = fSectorSumx2[sector];
     601
     602    // 1. Calculate the mean of the sums:
     603    Float_t ped        = sum/nevts;
     604
     605    // 2. Calculate the Variance of the sums:
     606    Float_t var = (sum2-sum*sum/nevts)/(nevts-1.);
     607
     608    // 3. Calculate the amplitude of the 150MHz "AB" noise
     609    Float_t abOffs = (fSectorSumAB0[sector] - fSectorSumAB1[sector]) / nevts;
     610
     611    // 4. Scale the mean, variance and AB-noise to the number of slices:
     612    ped    /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
     613    var    /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
     614    abOffs /= fExtractor ? fExtractor->GetNumHiGainSamples() : fExtractWinSize;
     615
     616    // 5. Scale the mean, variance and AB-noise to the number of pixels:
     617    ped    /= nspix;
     618    var    /= nspix;
     619    abOffs /= nspix;
     620
     621    // 6. Calculate the RMS from the Variance:
     622    const Float_t rms = var<0 ? 0 : TMath::Sqrt(var);
     623
     624    fPedestalsOut->GetAverageSector(sector).Set(ped, rms, abOffs, nevts);
     625}
    636626
    637627void MExtractPedestal::Print(Option_t *o) const
    638628{
    639 
    640629    *fLog << GetDescriptor() << ":" << endl;
    641     *fLog << "Name of incoming MPedestalCam Container:  " << fNamePedestalCamIn.Data() << endl;
    642     *fLog << "Name of outgoing MPedestalCam Container:  " << fNamePedestalCamOut.Data() << endl;
    643     *fLog << "Number events for pedestal calculation:  " << fNumEventsDump << endl;
    644     *fLog << "Number events for av. areas calculation: " << fNumAreasDump << endl;
    645     *fLog << "Number events for av. sector calculation: " << fNumSectorsDump << endl;
    646     *fLog << "Pedestal Update is                        " << (fPedestalUpdate?"on":"off") << endl;
     630    *fLog << "Name of input  MPedestalCam:  " << (fPedestalsIn?fPedestalsIn->GetName():fNamePedestalCamIn.Data()) << " (" << fPedestalsIn << ")" << endl;
     631    *fLog << "Name of output MPedestalCam:  " << (fPedestalsOut?fPedestalsOut->GetName():fNamePedestalCamOut.Data()) << " (" << fPedestalsOut << ")" << endl;
     632    *fLog << "Num evts for pedestal   calc: " << fNumEventsDump << endl;
     633    *fLog << "Num evts for avg.areas  calc: " << fNumAreasDump << endl;
     634    *fLog << "Num evts for avg.sector calc: " << fNumSectorsDump << endl;
     635    *fLog << "Pedestal Update is            " << (fPedestalUpdate?"on":"off") << endl;
    647636    *fLog << "ExtractWindow from slice " << fExtractWinFirst << " to " << fExtractWinLast << " incl." << endl;
    648 }
     637
     638    if (fExtractor)
     639        *fLog << "Extractor used:               " << fExtractor->ClassName() << endl;
     640}
  • trunk/MagicSoft/Mars/mpedestal/MExtractPedestal.h

    r5474 r5558  
    1414#endif
    1515
    16 #ifndef ROOT_TString
    17 #include "TString.h"
    18 #endif
    19 
    2016class MGeomCam;
    2117class MPedestalCam;
     
    2420class MRawEvtHeader;
    2521class MExtractTimeAndCharge;
     22
    2623class MExtractPedestal : public MTask
    2724{
    2825private:
     26  static const TString  fgNamePedestalCam;  //! "MPedestalCam"
    2927
    30   static const TString  fgNamePedestalCam;  //! "MPedestalCam"
    31    
    3228protected:
    3329 
     
    4945  UInt_t  fNumEventsDump;       // Number of event after which MPedestalCam gets updated
    5046  UInt_t  fNumAreasDump;        // Number of events after which averaged areas gets updated
    51   UInt_t  fNumSectorsDump;      // Number of events after which averaged sectors gets updated 
     47  UInt_t  fNumSectorsDump;      // Number of events after which averaged sectors gets updated
    5248
    5349  Bool_t  fPedestalUpdate;     // Flag if the pedestal shall be updated after every fNumEventsDump
     
    7066  MArrayI fSectorValid;        // number of valid pixels within sector idx
    7167
    72   Int_t  PreProcess (MParList *pList);
    73   Bool_t ReInit     (MParList *pList);
     68  Int_t  PreProcess(MParList *pList);
     69  Int_t  PostProcess();
     70  Bool_t ReInit(MParList *pList);
    7471  Int_t  ReadEnv(const TEnv &env, TString prefix, Bool_t print);
    7572
    7673  virtual void ResetArrays();
    77   void CalcPixResults   ( const UInt_t nevts, const UInt_t pixid );
    78   void CalcAreaResults  ( const UInt_t nevts, const UInt_t napix, const UInt_t aidx );
    79   void CalcSectorResults( const UInt_t nevts, const UInt_t nspix, const UInt_t sector ); 
    80  
     74
     75  void CalcPixResults   (const UInt_t nevts, const UInt_t pixid);
     76  void CalcAreaResults  (const UInt_t nevts, const UInt_t napix, const UInt_t aidx);
     77  void CalcSectorResults(const UInt_t nevts, const UInt_t nspix, const UInt_t sector);
     78
    8179public:
    82 
    8380  MExtractPedestal(const char *name=NULL, const char *title=NULL);
    8481
     
    8784
    8885  Bool_t SetExtractWindow(UShort_t first, UShort_t size);
    89   void SetNumEventsDump  ( UInt_t dumpevents = 500  ) { fNumEventsDump  = dumpevents; }
    90   void SetNumAreasDump   ( UInt_t dumpevents = 500  ) { fNumAreasDump   = dumpevents; } 
    91   void SetNumSectorsDump ( UInt_t dumpevents = 500  ) { fNumSectorsDump = dumpevents; }
    92   void SetPedestalUpdate ( Bool_t pedupdate         ) { fPedestalUpdate = pedupdate;  }
    93   void SetNamePedestalCamIn( const char *name=fgNamePedestalCam.Data()) { fNamePedestalCamIn = name; }
    94   void SetNamePedestalCamOut( const char *name=fgNamePedestalCam.Data()) { fNamePedestalCamOut = name; } 
    95   void SetExtractor      ( MExtractTimeAndCharge *e ) { fExtractor      = e; }
    96  
     86
     87  void SetNumEventsDump  (UInt_t dumpevents=500)  { fNumEventsDump  = dumpevents; }
     88  void SetNumAreasDump   (UInt_t dumpevents=500)  { fNumAreasDump   = dumpevents; }
     89  void SetNumSectorsDump (UInt_t dumpevents=500)  { fNumSectorsDump = dumpevents; }
     90  void SetPedestalUpdate (Bool_t pedupdate=kTRUE) { fPedestalUpdate = pedupdate;  }
     91
     92  void SetNamePedestalCamIn(const char *name=fgNamePedestalCam.Data())  { fNamePedestalCamIn = name; }
     93  void SetNamePedestalCamOut(const char *name=fgNamePedestalCam.Data()) { fNamePedestalCamOut = name; }
     94
     95  void SetPedestalsIn(MPedestalCam *pedcam) { fPedestalsIn = pedcam; }
     96
     97  void SetExtractor(MExtractTimeAndCharge *e) { fExtractor = e; }
     98
    9799  ClassDef(MExtractPedestal, 0)   // Base class for pedestal extractors
    98100};
  • trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.cc

    r5547 r5558  
    134134/////////////////////////////////////////////////////////////////////////////
    135135#include "MPedCalcFromLoGain.h"
    136 #include "MExtractPedestal.h"
    137 
    138 #include "MExtractTimeAndCharge.h"
    139136
    140137#include "MParList.h"
     
    152149#include "MGeomPix.h"
    153150#include "MGeomCam.h"
     151
     152#include "MExtractPedestal.h"
     153#include "MExtractTimeAndCharge.h"
    154154
    155155ClassImp(MPedCalcFromLoGain);
     
    162162const UShort_t MPedCalcFromLoGain::fgExtractWinSize  =  6;
    163163const UShort_t MPedCalcFromLoGain::fgMaxSignalVar    = 40;
     164
    164165// --------------------------------------------------------------------------
    165166//
     
    173174// - AddToBranchList("fHiGainFadcSamples");
    174175// - SetCheckRange(fgCheckWinFirst, fgCheckWinLast, fgExtractWinFirst, fgExtractWinSize)
    175 // - Clear()
    176176//
    177177MPedCalcFromLoGain::MPedCalcFromLoGain(const char *name, const char *title)
    178178{
    179 
    180   fName  = name  ? name  : "MPedCalcFromLoGain";
    181   fTitle = title ? title : "Task to calculate pedestals from lo-gains";
    182  
    183   SetCheckRange(fgCheckWinFirst, fgCheckWinLast);
    184   SetExtractWindow(fgExtractWinFirst, fgExtractWinSize);
    185  
    186   SetMaxSignalVar(fgMaxSignalVar);
    187  
    188   Clear();
     179    fName  = name  ? name  : "MPedCalcFromLoGain";
     180    fTitle = title ? title : "Task to calculate pedestals from lo-gains";
     181
     182    SetCheckRange(fgCheckWinFirst, fgCheckWinLast);
     183    SetExtractWindow(fgExtractWinFirst, fgExtractWinSize);
     184
     185    SetMaxSignalVar(fgMaxSignalVar);
    189186}
    190187
    191188void MPedCalcFromLoGain::ResetArrays()
    192189{
    193 
    194   MExtractPedestal::ResetArrays();
    195  
    196   fNumEventsUsed.Reset();
    197   fTotalCounter.Reset();
     190    MExtractPedestal::ResetArrays();
     191
     192    fNumEventsUsed.Reset();
     193    fTotalCounter.Reset();
    198194}
    199195
     
    235231Bool_t MPedCalcFromLoGain::ReInit(MParList *pList)
    236232{
    237 
    238   const UShort_t hisamples = fRunHeader->GetNumSamplesHiGain();
    239   const UShort_t losamples = fRunHeader->GetNumSamplesLoGain();
    240 
    241   UShort_t lastavailable   = hisamples+losamples-1;
    242 
    243   if (fExtractor)
    244     lastavailable = losamples-1;
    245 
    246   // If the size is not yet set, set the size
    247   if (fSumx.GetSize()==0)
    248   {
    249     const Int_t npixels  = fPedestalsOut->GetSize();
    250     fNumEventsUsed.Set(npixels);
    251     fTotalCounter.Set(npixels);
    252   }
    253  
    254   if (fExtractWinLast > lastavailable) //changed to override check
    255    {
    256      const UShort_t diff = fExtractWinLast - lastavailable;
    257      *fLog << warn << GetDescriptor();
    258      *fLog << " - WARNING: Selected Extract Window ranges out of range...adjusting to last available slice ";
    259      *fLog << lastavailable << endl;
    260      
    261      fExtractWinLast -= diff;
    262      fExtractWinSize -= diff;
    263    }
    264  
    265    lastavailable = hisamples + losamples -1;
    266    
    267    if (fCheckWinLast > lastavailable) //changed to override check
    268      {
    269        *fLog << warn << GetDescriptor()
    270              << " - WARNING: Last Check Window slice out of range...adjusting to last available slice "
    271              << lastavailable << endl;
    272        
    273        fCheckWinLast = lastavailable;
    274      }
    275 
    276    if (fCheckWinLast < fExtractWinLast)
    277      {
    278        *fLog << warn << GetDescriptor()
    279              << " - WARNING: Last Check Window slice lower than extraction window "
    280              << " - push it up to last number of low-gain samples: ";
    281        fCheckWinLast = fExtractWinLast;
    282        *fLog << fCheckWinLast << endl;
    283      }
    284 
    285   MExtractPedestal::ReInit(pList);
    286  
    287   return kTRUE;
    288 }
    289 
    290 // ---------------------------------------------------------------------------------
    291 //
    292 // Returns the pointer to slice "slice".
    293 //
    294 UShort_t MPedCalcFromLoGain::GetSlice(MRawEvtPixelIter *pixel, UInt_t slice)
    295 {
    296     const UShort_t nh = (Int_t)fRunHeader->GetNumSamplesHiGain();
    297 
    298     Byte_t *ptr;
    299 
    300     if(slice<nh)
    301         ptr = pixel->GetHiGainSamples() + slice;
    302     else
    303         ptr = pixel->GetLoGainSamples() + slice - nh;
    304 
    305    return *ptr;
     233    const UShort_t hisamples = fRunHeader->GetNumSamplesHiGain();
     234    const UShort_t losamples = fRunHeader->GetNumSamplesLoGain();
     235 
     236    fSlices.Set(hisamples+losamples);
     237 
     238    UShort_t lastavailable   = hisamples+losamples-1;
     239 
     240    if (fExtractor)
     241        lastavailable = losamples-1;
     242 
     243    // If the size is not yet set, set the size
     244    if (fSumx.GetSize()==0)
     245    {
     246        const Int_t npixels  = fPedestalsOut->GetSize();
     247        fNumEventsUsed.Set(npixels);
     248        fTotalCounter.Set(npixels);
     249    }
     250 
     251    if (fExtractWinLast > lastavailable) //changed to override check
     252    {
     253        const UShort_t diff = fExtractWinLast - lastavailable;
     254        *fLog << warn << GetDescriptor();
     255        *fLog << " - WARNING: Selected Extract Window ranges out of range...adjusting to last available slice ";
     256        *fLog << lastavailable << endl;
     257 
     258        fExtractWinLast -= diff;
     259        fExtractWinSize -= diff;
     260    }
     261 
     262    lastavailable = hisamples + losamples -1;
     263 
     264    if (fCheckWinLast > lastavailable) //changed to override check
     265    {
     266        *fLog << warn << GetDescriptor();
     267        *fLog << " - WARNING: Last Check Window slice out of range...adjusting to last available slice ";
     268        *fLog << lastavailable << endl;
     269 
     270        fCheckWinLast = lastavailable;
     271    }
     272 
     273    return MExtractPedestal::ReInit(pList);
    306274}
    307275
     
    314282Int_t MPedCalcFromLoGain::Process()
    315283{
    316 
    317   MRawEvtPixelIter pixel(fRawEvt);
    318  
    319   while (pixel.Next())
    320     {
    321       const UInt_t idx = pixel.GetPixelId();
    322       const UInt_t aidx   = (*fGeom)[idx].GetAidx();
    323       const UInt_t sector = (*fGeom)[idx].GetSector();     
    324      
    325       UShort_t max = 0;
    326       UShort_t min = (UShort_t)-1;
    327      
    328       // Find the maximum and minimum signal per slice in the high gain window
    329       for (Int_t slice=fCheckWinFirst; slice<=fCheckWinLast; slice++)
     284    // This is the workaround to put hi- and lo-gains together
     285    const Int_t nhigain = fRunHeader->GetNumSamplesHiGain();
     286    const Int_t nlogain = fRunHeader->GetNumSamplesLoGain();
     287
     288    Byte_t *slices = fSlices.GetArray();
     289
     290    // Real Process
     291    MRawEvtPixelIter pixel(fRawEvt);
     292
     293    while (pixel.Next())
     294    {
     295        // This is the fast workaround to put hi- and lo-gains together
     296        memcpy(slices,         pixel.GetHiGainSamples(), nhigain);
     297        memcpy(slices+nhigain, pixel.GetLoGainSamples(), nlogain);
     298
     299        // Start 'real' work
     300        const UInt_t idx = pixel.GetPixelId();
     301
     302        const UInt_t aidx   = (*fGeom)[idx].GetAidx();
     303        const UInt_t sector = (*fGeom)[idx].GetSector();
     304
     305        UShort_t max = 0;
     306        UShort_t min = (UShort_t)-1;
     307
     308        // Find the maximum and minimum signal per slice in the high gain window
     309        for (Byte_t *slice=slices+fCheckWinFirst; slice<=slices+fCheckWinLast; slice++)
    330310        {
    331           const UShort_t svalue = GetSlice(&pixel,slice);
    332           if (svalue > max)
    333             max = svalue;
    334           if (svalue < min)
    335             min = svalue;
     311            if (*slice > max)
     312                max = *slice;
     313            if (*slice < min)
     314                min = *slice;
    336315        }
    337316     
    338       // If the maximum in the high gain window is smaller than
    339       if (max-min>=fMaxSignalVar || max>=250)
    340         continue;
    341      
    342       Float_t sum  = 0.;
    343       UInt_t  sumi = 0;
    344      
    345       //extract pedestal
    346       if (fExtractor)
    347         CalcExtractor( &pixel, sum, (*fPedestalsIn)[idx] );
    348       else
     317        // If the maximum in the high gain window is smaller than
     318        if (max-min>=fMaxSignalVar || max>=250)
     319            continue;
     320
     321        Float_t sum  = 0.;
     322        UInt_t  sumi = 0;
     323
     324        //extract pedestal
     325        if (fExtractor)
     326            CalcExtractor(pixel, sum, (*fPedestalsIn)[idx]);
     327        else
    349328        {
    350           for(Int_t slice=fExtractWinFirst; slice<=fExtractWinLast; slice++)
    351             sumi += GetSlice(&pixel,slice);
    352           sum = (Float_t)sumi;
     329            for(Byte_t *slice=slices+fExtractWinFirst; slice<=slices+fExtractWinLast; slice++)
     330                sumi += *slice;
     331            sum = (Float_t)sumi;
    353332        }
    354      
    355       const Float_t sqrsum = sum*sum;
    356      
    357       fSumx[idx]           += sum;
    358       fSumx2[idx]          += sqrsum;
    359       fAreaSumx[aidx]      += sum;
    360       fAreaSumx2[aidx]     += sqrsum;
    361       fSectorSumx[sector]  += sum;     
    362       fSectorSumx2[sector] += sqrsum;     
    363      
    364       fNumEventsUsed[idx]   ++;
    365       fAreaFilled   [aidx]  ++;
    366       fSectorFilled [sector]++;
    367      
    368       if (!fExtractor)
     333
     334        const Float_t sqrsum = sum*sum;
     335
     336        fSumx[idx]           += sum;
     337        fSumx2[idx]          += sqrsum;
     338        fAreaSumx[aidx]      += sum;
     339        fAreaSumx2[aidx]     += sqrsum;
     340        fSectorSumx[sector]  += sum;
     341        fSectorSumx2[sector] += sqrsum;
     342
     343        fNumEventsUsed[idx]   ++;
     344        fAreaFilled   [aidx]  ++;
     345        fSectorFilled [sector]++;
     346
     347        if (!fExtractor)
    369348        {
    370           //
    371           // Calculate the amplitude of the 150MHz "AB" noise
    372           //
    373           const UShort_t abFlag = (fExtractWinFirst + pixel.HasABFlag()) & 0x1;
    374          
    375           for (UShort_t islice=fExtractWinFirst; islice<=fExtractWinLast; islice+=2)
     349            //
     350            // Calculate the amplitude of the 150MHz "AB" noise
     351            //
     352            const UShort_t abFlag = (fExtractWinFirst + pixel.HasABFlag()) & 0x1;
     353
     354            for (Byte_t *slice=slices+fExtractWinFirst; slice<=slices+fExtractWinLast; slice+=2)
    376355            {
    377               const UShort_t sliceAB0 = islice + abFlag;
    378               const UShort_t sliceAB1 = islice - abFlag + 1;
    379               const UShort_t ab0 = GetSlice(&pixel, sliceAB0);
    380               const UShort_t ab1 = GetSlice(&pixel, sliceAB1);
    381               fSumAB0[idx]        += ab0;
    382               fSumAB1[idx]        += ab1;
    383               fAreaSumAB0[aidx]   += ab0;
    384               fAreaSumAB1[aidx]   += ab1;     
    385               fSectorSumAB0[aidx] += ab0;
    386               fSectorSumAB1[aidx] += ab1;     
     356                const UShort_t ab0 = *(slice + abFlag);
     357                const UShort_t ab1 = *(slice - abFlag + 1);
     358
     359                fSumAB0[idx]        += ab0;
     360                fSumAB1[idx]        += ab1;
     361                fAreaSumAB0[aidx]   += ab0;
     362                fAreaSumAB1[aidx]   += ab1;
     363                fSectorSumAB0[aidx] += ab0;
     364                fSectorSumAB1[aidx] += ab1;
    387365            }
    388366        }
    389      
    390       if (!fPedestalUpdate || (UInt_t)fNumEventsUsed[idx]<fNumEventsDump)
    391         continue;
    392      
    393       CalcPixResults(fNumEventsDump, idx);
    394       fTotalCounter[idx]++;
    395      
    396       fNumEventsUsed[idx]=0;
    397       fSumx[idx]=0;
    398       fSumx2[idx]=0;
    399       fSumAB0[idx]=0;
    400       fSumAB1[idx]=0;
    401     }
    402  
    403   if (!(GetNumExecutions() % fNumAreasDump))
    404     {
    405       //
    406       // Loop over the (two) area indices to get the averaged pedestal per aidx
    407       //
    408       for (UInt_t aidx=0; aidx<fAreaValid.GetSize(); aidx++)
     367
     368        if (!fPedestalUpdate || (UInt_t)fNumEventsUsed[idx]<fNumEventsDump)
     369            continue;
     370
     371        CalcPixResults(fNumEventsDump, idx);
     372        fTotalCounter[idx]++;
     373
     374        fNumEventsUsed[idx]=0;
     375        fSumx[idx]=0;
     376        fSumx2[idx]=0;
     377        fSumAB0[idx]=0;
     378        fSumAB1[idx]=0;
     379    }
     380
     381    if (!(GetNumExecutions() % fNumAreasDump))
     382        CalcAreaResult();
     383
     384    if (!(GetNumExecutions() % fNumSectorsDump))
     385        CalcSectorResult();
     386
     387    if (fPedestalUpdate)
     388        fPedestalsOut->SetReadyToSave();
     389
     390  return kTRUE;
     391}
     392
     393// --------------------------------------------------------------------------
     394//
     395// Loop over the sector indices to get the averaged pedestal per sector
     396//
     397void MPedCalcFromLoGain::CalcSectorResult()
     398{
     399    for (UInt_t sector=0; sector<fSectorFilled.GetSize(); sector++)
     400        if (fSectorValid[sector]>0)
     401            CalcSectorResults(fSectorFilled[sector], fSectorValid[sector], sector);
     402}
     403
     404// --------------------------------------------------------------------------
     405//
     406// Loop over the (two) area indices to get the averaged pedestal per aidx
     407//
     408void MPedCalcFromLoGain::CalcAreaResult()
     409{
     410    for (UInt_t aidx=0; aidx<fAreaFilled.GetSize(); aidx++)
     411        if (fAreaValid[aidx]>0)
     412            CalcAreaResults(fAreaFilled[aidx], fAreaValid[aidx], aidx);
     413
     414}
     415
     416void MPedCalcFromLoGain::CalcExtractor(const MRawEvtPixelIter &pixel, Float_t &sum, MPedestalPix &ped)
     417{
     418    Byte_t  sat  = 0;
     419    Byte_t *logain = pixel.GetLoGainSamples() + fExtractWinFirst;
     420
     421    const Bool_t logainabflag = (pixel.HasABFlag() + pixel.GetNumHiGainSamples()) & 0x1;
     422
     423    Float_t dummy;
     424    fExtractor->FindTimeAndChargeHiGain(logain,logain,sum,dummy,dummy,dummy,sat,ped,logainabflag);
     425}
     426
     427
     428// --------------------------------------------------------------------------
     429//
     430// Compute signal mean and rms in the whole run and store it in MPedestalCam
     431//
     432Int_t MPedCalcFromLoGain::PostProcess()
     433{
     434    // Compute pedestals and rms from the whole run
     435    if (fPedestalUpdate)
     436        return kTRUE;
     437
     438    *fLog << flush << inf << "Calculating Pedestals..." << flush;
     439
     440    const Int_t npix = fGeom->GetNumPixels();
     441    for (Int_t idx=0; idx<npix; idx++)
     442    {
     443        const ULong_t n = fNumEventsUsed[idx];
     444        if (n>1)
    409445        {
    410           const Int_t   napix = fAreaValid.At(aidx);
    411           if (napix == 0)
    412             continue;
    413          
    414           CalcAreaResults(fSectorFilled[aidx],napix,aidx);
     446            CalcPixResults(n, idx);
     447            fTotalCounter[idx]++;
    415448        }
    416449    }
    417  
    418   if (!(GetNumExecutions() % fNumSectorsDump))
    419     {
    420       //
    421       // Loop over the (two) sector indices to get the averaged pedestal per sector
    422       //
    423       for (UInt_t sector=0; sector<fSectorValid.GetSize(); sector++)
    424         {
    425           const Int_t   nspix = fSectorValid.At(sector);
    426           if (nspix == 0)
    427             continue;
    428          
    429           CalcSectorResults(fSectorFilled[sector],nspix,sector);
    430         }
    431     }
    432  
    433   if (fPedestalUpdate)
     450
     451    CalcAreaResult();
     452    CalcSectorResult();
     453
    434454    fPedestalsOut->SetReadyToSave();
    435  
    436   return kTRUE;
    437 }
    438 
    439 void MPedCalcFromLoGain::CalcExtractor( MRawEvtPixelIter *pixel, Float_t &sum, MPedestalPix &ped)
    440 {
    441  
    442   Byte_t  sat  = 0;
    443   Byte_t *logain = pixel->GetLoGainSamples() + fExtractWinFirst;
    444   const Bool_t logainabflag = (pixel->HasABFlag() + pixel->GetNumHiGainSamples()) & 0x1;
    445   Float_t dummy;
    446   fExtractor->FindTimeAndChargeHiGain(logain,logain,sum,dummy,dummy,dummy,sat,ped,logainabflag);
    447 }
    448 
    449 
    450 // --------------------------------------------------------------------------
    451 //
    452 // Compute signal mean and rms in the whole run and store it in MPedestalCam
    453 //
    454 Int_t MPedCalcFromLoGain::PostProcess()
    455 {
    456 
    457   // Compute pedestals and rms from the whole run
    458   if (fPedestalUpdate)
    459     return kTRUE;
    460  
    461   *fLog << flush << inf << "Calculating Pedestals..." << flush;
    462  
    463   const Int_t npix = fGeom->GetNumPixels();
    464   for (Int_t idx=0; idx<npix; idx++)
    465     {
    466       const ULong_t n = fNumEventsUsed[idx];
    467       if (n>1)
    468         {
    469           CalcPixResults(n, idx);
    470           fTotalCounter[idx]++;
    471         }
    472     }
    473 
    474   //
    475   // Loop over the (two) area indices to get the averaged pedestal per aidx
    476   //
    477   for (UInt_t aidx=0; aidx<fAreaFilled.GetSize(); aidx++)
    478     {
    479       const Int_t   napix = fAreaValid.At(aidx);
    480       if (napix == 0)
    481         continue;
    482 
    483       CalcAreaResults(fAreaFilled[aidx],napix,aidx);
    484     }
    485  
    486   //
    487   // Loop over the (six) sector indices to get the averaged pedestal per sector
    488   //
    489   for (UInt_t sector=0; sector<fSectorFilled.GetSize(); sector++)
    490     {
    491       const Int_t   nspix = fSectorValid.At(sector);
    492       if (nspix == 0)
    493         continue;
    494      
    495       CalcSectorResults(fSectorFilled[sector],nspix,sector);
    496     }
    497 
    498   fPedestalsOut->SetReadyToSave();
    499 
    500   return kTRUE;
     455
     456    return MExtractPedestal::PostProcess();
    501457}
    502458
     
    511467Int_t MPedCalcFromLoGain::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
    512468{
    513 
    514   Bool_t rc=kFALSE;
    515  
    516   // Find resources for CheckWindow
    517   Int_t fs = fCheckWinFirst;
    518   Int_t ls = fCheckWinLast;
    519   if (IsEnvDefined(env, prefix, "CheckWinFirst", print))
    520     {
    521       fs = GetEnvValue(env, prefix, "CheckWinFirst", fs);
    522       rc = kTRUE;
    523     }
    524   if (IsEnvDefined(env, prefix, "CheckWinLast", print))
    525     {
    526       ls = GetEnvValue(env, prefix, "CheckWinLast", ls);
    527       rc = kTRUE;
    528     }
    529  
    530   SetCheckRange(fs,ls);
    531  
    532   // find resource for maximum signal variation
    533   if (IsEnvDefined(env, prefix, "MaxSignalVar", print))
    534     {
    535       SetMaxSignalVar(GetEnvValue(env, prefix, "MaxSignalVar", fMaxSignalVar));
    536       rc = kTRUE;
    537     }
    538  
    539   return MExtractPedestal::ReadEnv(env,prefix,print) ? kTRUE : rc;
     469    Bool_t rc=kFALSE;
     470
     471    // Find resources for CheckWindow
     472    Int_t fs = fCheckWinFirst;
     473    Int_t ls = fCheckWinLast;
     474    if (IsEnvDefined(env, prefix, "CheckWinFirst", print))
     475    {
     476        fs = GetEnvValue(env, prefix, "CheckWinFirst", fs);
     477        rc = kTRUE;
     478    }
     479    if (IsEnvDefined(env, prefix, "CheckWinLast", print))
     480    {
     481        ls = GetEnvValue(env, prefix, "CheckWinLast", ls);
     482        rc = kTRUE;
     483    }
     484
     485    SetCheckRange(fs,ls);
     486
     487    // find resource for maximum signal variation
     488    if (IsEnvDefined(env, prefix, "MaxSignalVar", print))
     489    {
     490        SetMaxSignalVar(GetEnvValue(env, prefix, "MaxSignalVar", fMaxSignalVar));
     491        rc = kTRUE;
     492    }
     493
     494    return MExtractPedestal::ReadEnv(env,prefix,print) ? kTRUE : rc;
    540495}
    541496
    542497void MPedCalcFromLoGain::Print(Option_t *o) const
    543498{
    544 
    545   MExtractPedestal::Print(o);
    546 
    547   *fLog << "CheckWindow   from slice  " << fCheckWinFirst   << " to " << fCheckWinLast << " incl." << endl;
    548   *fLog << "Max. allowed signal variation:             " << fMaxSignalVar << endl;
    549   *fLog << endl;
    550 }
     499    MExtractPedestal::Print(o);
     500
     501    *fLog << "CheckWindow   from slice  " << fCheckWinFirst   << " to " << fCheckWinLast << " incl." << endl;
     502    *fLog << "Max. allowed signal variation:             " << fMaxSignalVar << endl;
     503    *fLog << endl;
     504}
  • trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.h

    r5498 r5558  
    66#endif
    77
    8 #ifndef ROOT_TArrayD
    9 #include <TArrayD.h>
    10 #endif
    11 
    128#ifndef ROOT_TArrayI
    139#include <TArrayI.h>
    1410#endif
    1511
     12#ifndef ROOT_MArrayB
     13#include "MArrayB.h"
     14#endif
     15
    1616class MRawEvtPixelIter;
    1717class MPedestalPix;
     18
    1819class MPedCalcFromLoGain : public MExtractPedestal
    1920{
     21private:
     22    static const UShort_t fgCheckWinFirst;    // First FADC slice to check for signal (currently set to: 0)
     23    static const UShort_t fgCheckWinLast;     // Last FADC slice to check for signal  (currently set to: 29)
     24    static const UShort_t fgMaxSignalVar;     // The maximum difference between the highest and lowest slice
     25    static const UShort_t fgExtractWinFirst;  // First FADC slice to use for pedestal calculation (currently set to: 15)
     26    static const UShort_t fgExtractWinSize;   // number of successive slices used to calculate pedestal (currently set to: 6)
    2027
    21   static const UShort_t fgCheckWinFirst;    // First FADC slice to check for signal (currently set to: 0)
    22   static const UShort_t fgCheckWinLast;     // Last FADC slice to check for signal  (currently set to: 29)
    23   static const UShort_t fgMaxSignalVar;     // The maximum difference between the highest and lowest slice
    24   static const UShort_t fgExtractWinFirst;  // First FADC slice to use for pedestal calculation (currently set to: 15)
    25   static const UShort_t fgExtractWinSize;   // number of successive slices used to calculate pedestal (currently set to: 6)
    26    
    27   UShort_t fMaxSignalVar;
    28   UShort_t fCheckWinFirst;
    29   UShort_t fCheckWinLast;
     28    UShort_t fMaxSignalVar;
     29    UShort_t fCheckWinFirst;
     30    UShort_t fCheckWinLast;
    3031
    31   TArrayI fNumEventsUsed;      // Number of events used for pedestal calc for each pixel
    32   TArrayI fTotalCounter;       // Counter for dumping values to Pedestal Container
     32    TArrayI fNumEventsUsed;      //! Number of events used for pedestal calc for each pixel
     33    TArrayI fTotalCounter;       //! Counter for dumping values to Pedestal Container
    3334
    34   Bool_t ReInit     (MParList *pList);
    35   Int_t  Process    ();
    36   Int_t  PostProcess();
    37  
    38   Int_t  ReadEnv(const TEnv &env, TString prefix, Bool_t print);
    39  
    40   //Helper function to extract slice values by slice number
    41   UShort_t GetSlice(MRawEvtPixelIter *pixel, UInt_t slice);
    42   void CalcExtractor   ( MRawEvtPixelIter *pixel, Float_t &sum, MPedestalPix &ped);
    43   void ResetArrays();
     35    MArrayB fSlices;             //! workaround to put hi- and lo-gain slices together
     36
     37    Bool_t ReInit(MParList *pList);
     38    Int_t  Process();
     39    Int_t  PostProcess();
     40
     41    Int_t  ReadEnv(const TEnv &env, TString prefix, Bool_t print);
     42
     43    //Helper function to extract slice values by slice number
     44    void CalcExtractor(const MRawEvtPixelIter &pixel, Float_t &sum, MPedestalPix &ped);
     45    void ResetArrays();
     46
     47    void CalcSectorResult();
     48    void CalcAreaResult();
    4449
    4550public:
    4651
    47   MPedCalcFromLoGain(const char *name=NULL, const char *title=NULL);
    48  
    49   void Print(Option_t *o="") const;
    50  
    51   // Setters
    52   Bool_t SetCheckRange(UShort_t checkfirst=fgCheckWinFirst, UShort_t checklast=fgCheckWinLast);
    53   void SetMaxSignalVar(UShort_t maxvar=40)       { fMaxSignalVar = maxvar;    }
     52    MPedCalcFromLoGain(const char *name=NULL, const char *title=NULL);
    5453
    55   // Getters
    56   TArrayI *GetNumEventsUsed() { return &fNumEventsUsed; }
    57  
    58   ClassDef(MPedCalcFromLoGain, 1)   // Task to calculate pedestals from data runs
     54    void Print(Option_t *o="") const;
     55
     56    // Setters
     57    Bool_t SetCheckRange(UShort_t checkfirst=fgCheckWinFirst, UShort_t checklast=fgCheckWinLast);
     58    void SetMaxSignalVar(UShort_t maxvar=40)       { fMaxSignalVar = maxvar;    }
     59
     60    // Getters
     61    TArrayI *GetNumEventsUsed() { return &fNumEventsUsed; }
     62
     63    ClassDef(MPedCalcFromLoGain, 1)   // Task to calculate pedestals from data runs
    5964};
    6065
  • trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.cc

    r5550 r5558  
    119119/////////////////////////////////////////////////////////////////////////////
    120120#include "MPedCalcPedRun.h"
    121 #include "MExtractPedestal.h"
    122 
    123 #include "MExtractTimeAndCharge.h"
    124121
    125122#include "MParList.h"
     
    139136#include "MGeomCam.h"
    140137
    141 #include "MStatusDisplay.h"
     138#include "MExtractPedestal.h"
     139#include "MExtractTimeAndCharge.h"
    142140
    143141ClassImp(MPedCalcPedRun);
     
    148146const UShort_t MPedCalcPedRun::fgExtractWinSize        = 6;
    149147const UInt_t   MPedCalcPedRun::gkFirstRunWithFinalBits = 45589;
     148
    150149// --------------------------------------------------------------------------
    151150//
     
    161160// - AddToBranchList("fHiGainFadcSamples");
    162161// - SetRange(fgHiGainFirst, fgHiGainLast, fgLoGainFirst, fgLoGainLast)
    163 // - Clear()
    164162//
    165163MPedCalcPedRun::MPedCalcPedRun(const char *name, const char *title)
    166 {
    167   fName  = name  ? name  : "MPedCalcPedRun";
    168   fTitle = title ? title : "Task to calculate pedestals from pedestal runs raw data";
    169 
    170   AddToBranchList("fHiGainPixId");
    171   AddToBranchList("fHiGainFadcSamples");
    172  
    173   SetExtractWindow(fgExtractWinFirst, fgExtractWinSize);
    174 
    175   Clear();
    176  
    177   fFirstRun   = kTRUE;
    178   fSkip       = kFALSE;
    179   fOverlap    = 0;
    180   fUsedEvents = 0;
    181 }
    182 
    183 
    184 
    185 void MPedCalcPedRun::ResetArrays()
    186 {
    187 
    188   MExtractPedestal::ResetArrays();
    189   fUsedEvents    = 0;
     164    : fOverlap(0)
     165{
     166    fName  = name  ? name  : "MPedCalcPedRun";
     167    fTitle = title ? title : "Task to calculate pedestals from pedestal runs raw data";
     168
     169    SetExtractWindow(fgExtractWinFirst, fgExtractWinSize);
    190170}
    191171
     
    197177// same number of slices.
    198178//
    199 // The run type is checked for "Pedestal" (run type: 1)
     179void MPedCalcPedRun::CheckExtractionWindow()
     180{
     181    const UShort_t lastavailable = fRunHeader->GetNumSamplesHiGain()-1;
     182
     183    if (fExtractWinLast <= lastavailable)
     184        return;
     185
     186    const UShort_t diff = fExtractWinLast - lastavailable;
     187
     188    *fLog << warn << endl;
     189    *fLog << "Selected ExtractWindow [" << fExtractWinFirst << "," << fExtractWinLast;
     190    *fLog << "] ranges out of range [0," << lastavailable << "]." << endl;
     191    *fLog << "Using " << diff << " samples from the 'Lo-Gain' slices for the pedestal extraction" << endl;
     192
     193    fExtractWinLast -= diff;
     194    fOverlap         = diff;
     195}
     196
     197// --------------------------------------------------------------------------
     198//
     199// Set fIsFirstPedRun=kTRUE
     200//
     201Int_t MPedCalcPedRun::PreProcess(MParList *pList)
     202{
     203    fUsedEvents = 0;
     204    fIsFirstPedRun = kTRUE;
     205    return MExtractPedestal::PreProcess(pList);
     206}
     207
     208// --------------------------------------------------------------------------
     209//
     210// The run type is checked for "kRTPedestal"
    200211// and the variable fSkip is set in that case
    201212//
    202213Bool_t MPedCalcPedRun::ReInit(MParList *pList)
    203214{
    204 
    205   MExtractPedestal::ReInit(pList);
    206 
    207   const UShort_t lastavailable = fRunHeader->GetNumSamplesHiGain()-1;
    208  
    209   if (fExtractWinLast > lastavailable)
     215    if (!MExtractPedestal::ReInit(pList))
     216        return kFALSE;
     217
     218    CheckExtractionWindow();
     219
     220    // If this is the first ped run, the next run (call to ReInit)
     221    // is not the first anymore
     222    if (fRunHeader->GetRunType()==MRawRunHeader::kRTPedestal)
    210223    {
    211       const UShort_t diff = fExtractWinLast - lastavailable;
    212       *fLog << endl;
    213       *fLog << warn << GetDescriptor()
    214             << Form("%s%2i%s%2i%s%2i%s",": Selected ExtractWindow [",
    215                     fExtractWinFirst,",",fExtractWinLast,
    216                     "] ranges out of available limits: [0,",lastavailable,"].") << endl;
    217       *fLog << warn << GetDescriptor()
    218             << Form("%s%2i%s",": Will use ",diff," samples from the 'Low-Gain' slices for the pedestal extraction")
    219             << endl;
    220       fExtractWinLast -= diff;
    221       fOverlap         = diff;
     224        fIsFirstPedRun = kFALSE;
     225        return kTRUE;
    222226    }
    223227
    224   const UShort_t runtype = fRunHeader->GetRunType();
    225 
    226   switch (runtype)
    227     {
    228     case 1:
    229       fFirstRun = kFALSE;
    230       fSkip     = kFALSE;
    231       return kTRUE;
    232       break;
    233     default:
    234       fSkip     = kTRUE;
    235       if (!fFirstRun)
    236         {
    237           *fLog << endl;
    238           *fLog << inf << GetDescriptor() << " : Finalize pedestal calculations..." << flush;
    239           Finalize();
    240           Reset();
    241         }
    242       return kTRUE;       
    243       break;
    244     }
    245  
    246   Print();
    247 
    248   return kTRUE;
     228    // If this is the first call to ReInit (before reading the first file)
     229    // nothing should be done
     230    if (fIsFirstPedRun)
     231        return kTRUE;
     232
     233    // In any other case some kind of finaliazation must be done
     234    *fLog << inf << "Finalizing pedestal calculations..." << flush;
     235
     236    if (!Finalize())
     237        return kFALSE;
     238
     239    Reset();
     240
     241    return kTRUE;
    249242}
    250243
     
    260253Int_t MPedCalcPedRun::Process()
    261254{
     255    if (!IsPedBitSet())
     256    {
     257        *fLog << err << GetDescriptor() << " - ERROR: IsPedBitSet() returned kFALSE... abort." << endl;
     258        return kERROR;
     259    }
     260
     261    if (fRunHeader->GetRunType() != MRawRunHeader::kRTPedestal)
     262        return kTRUE;
    262263
    263264  fUsedEvents++;
     
    266267 
    267268  while (pixel.Next())
    268     {
     269  {
    269270      const UInt_t idx    = pixel.GetPixelId();
    270271      const UInt_t aidx   = (*fGeom)[idx].GetAidx();
     
    276277     
    277278      if (fExtractor)
    278         CalcExtractor( &pixel, sum, (*fPedestalsIn)[idx]);
     279          CalcExtractor(pixel, sum, (*fPedestalsIn)[idx]);
    279280      else
    280         CalcSums( &pixel, sum, ab0, ab1);
    281 
    282       fSumx[idx]          += sum;
    283       fAreaSumx[aidx]     += sum;
    284       fSectorSumx[sector] += sum;     
    285 
    286       const Float_t sqrsum  = sum*sum;
    287       fSumx2[idx]          += sqrsum;
    288       fAreaSumx2[aidx]     += sqrsum;
    289       fSectorSumx2[sector] += sqrsum;     
    290 
    291       fSumAB0[idx]        += ab0;
    292       fSumAB1[idx]        += ab1;
    293 
    294       fAreaSumAB0[aidx]   += ab0;
    295       fAreaSumAB1[aidx]   += ab1;     
     281          CalcSums(pixel, sum, ab0, ab1);
     282
     283      fSumx[idx]           += sum;
     284      fAreaSumx[aidx]      += sum;
     285      fSectorSumx[sector]  += sum;
     286
     287      const Float_t sqrsum   = sum*sum;
     288      fSumx2[idx]           += sqrsum;
     289      fAreaSumx2[aidx]      += sqrsum;
     290      fSectorSumx2[sector]  += sqrsum;
     291
     292      fSumAB0[idx]         += ab0;
     293      fSumAB1[idx]         += ab1;
     294
     295      fAreaSumAB0[aidx]    += ab0;
     296      fAreaSumAB1[aidx]    += ab1;
    296297     
    297       fSectorSumAB0[aidx] += ab0;
    298       fSectorSumAB1[aidx] += ab1;     
    299     }
     298      fSectorSumAB0[aidx]  += ab0;
     299      fSectorSumAB1[aidx]  += ab1;
     300  }
    300301 
    301302  fPedestalsOut->SetReadyToSave();
     
    305306
    306307
    307 void MPedCalcPedRun::CalcExtractor( MRawEvtPixelIter *pixel, Float_t &sum, MPedestalPix &ped)
    308 {
    309  
    310   Byte_t  sat = 0;
    311   Byte_t *first  = pixel->GetHiGainSamples() + fExtractWinFirst;
    312   Byte_t *logain = pixel->GetLoGainSamples();
    313   const Bool_t abflag  = pixel->HasABFlag();
    314   Float_t dummy;
    315   fExtractor->FindTimeAndChargeHiGain(first,logain,sum,dummy,dummy,dummy,sat,ped,abflag);
    316 
    317 }
    318 
    319 
    320 void MPedCalcPedRun::CalcSums( MRawEvtPixelIter *pixel, Float_t &sum, UInt_t &ab0, UInt_t &ab1)
    321 {
    322 
    323   Byte_t *higain = pixel->GetHiGainSamples() + fExtractWinFirst;
     308void MPedCalcPedRun::CalcExtractor(const MRawEvtPixelIter &pixel, Float_t &sum, MPedestalPix &ped)
     309{
     310    const Bool_t abflag = pixel.HasABFlag();
     311
     312    Byte_t *first  = pixel.GetHiGainSamples() + fExtractWinFirst;
     313    Byte_t *logain = pixel.GetLoGainSamples();
     314
     315    Byte_t  sat = 0;
     316
     317    Float_t dummy;
     318    fExtractor->FindTimeAndChargeHiGain(first,logain,sum,dummy,dummy,dummy,sat,ped,abflag);
     319}
     320
     321void MPedCalcPedRun::CalcSums(const MRawEvtPixelIter &pixel, Float_t &sum, UInt_t &ab0, UInt_t &ab1)
     322{
     323
     324  Byte_t *higain = pixel.GetHiGainSamples() + fExtractWinFirst;
    324325  Byte_t *ptr = higain;
    325326  Byte_t *end = ptr + fExtractWinSize;
    326327
    327   const Bool_t abflag = pixel->HasABFlag();
     328  const Bool_t abflag = pixel.HasABFlag();
    328329
    329330  Int_t sumi = 0;
     
    331332  Int_t  cnt = 0;
    332333  do
    333     {
     334  {
    334335      sumi += *ptr;
    335       if (pixel->IsABFlagValid())
    336         {
     336      if (!pixel.IsABFlagValid())
     337          continue;
     338
     339      const Int_t abFlag = (fExtractWinFirst + abflag + cnt) & 0x1;
     340      if (abFlag)
     341          ab1 += *ptr;
     342      else
     343          ab0 += *ptr;
     344
     345      cnt++;
     346  } while (++ptr != end);
     347
     348  if (fOverlap != 0)
     349  {
     350      ptr = pixel.GetLoGainSamples();
     351      end = ptr + fOverlap;
     352
     353      do
     354      {
     355          sumi += *ptr;
     356          if (!pixel.IsABFlagValid())
     357              continue;
     358
    337359          const Int_t abFlag = (fExtractWinFirst + abflag + cnt) & 0x1;
    338360          if (abFlag)
    339             ab1 += *ptr;
     361              ab1 += *ptr;
    340362          else
    341             ab0 += *ptr;
    342          
     363              ab0 += *ptr;
     364
    343365          cnt++;
    344         }
    345     }
    346   while (++ptr != end);
    347  
    348   if (fOverlap != 0)
    349     {
    350       ptr = pixel->GetLoGainSamples();
    351       end = ptr + fOverlap;
    352      
    353       do
    354         {
    355           sumi += *ptr;
    356           if (pixel->IsABFlagValid())
    357             {
    358               const Int_t abFlag = (fExtractWinFirst + abflag + cnt) & 0x1;
    359               if (abFlag)
    360                 ab1 += *ptr;
    361               else
    362                 ab0 += *ptr;
    363              
    364               cnt++;
    365             }
    366         }
    367       while (++ptr != end);
    368     }
     366      } while (++ptr != end);
     367  }
    369368
    370369  sum = (Float_t)sumi;
    371 
    372370}
    373371
     
    382380// Compute signal mean and rms in the whole run and store it in MPedestalCam
    383381//
    384 void MPedCalcPedRun::Finalize()
    385 {
    386 
    387   if (fUsedEvents == 0)
    388     return;
    389 
    390   MRawEvtPixelIter pixel(fRawEvt);
    391  
    392   while (pixel.Next())
    393     {
    394       const Int_t  pixid  = pixel.GetPixelId();
    395       CalcPixResults(fUsedEvents,pixid);
    396     }
    397 
    398   //
    399   // Loop over the (two) area indices to get the averaged pedestal per aidx
    400   //
    401   for (UInt_t aidx=0; aidx<fAreaValid.GetSize(); aidx++)
    402     {
    403       const Int_t   napix = fAreaValid.At(aidx);
    404       if (napix == 0)
    405         continue;
    406 
    407       CalcAreaResults(fUsedEvents,napix,aidx);
    408     }
    409  
    410   //
    411   // Loop over the (six) sector indices to get the averaged pedestal per sector
    412   //
    413   for (UInt_t sector=0; sector<fSectorValid.GetSize(); sector++)
    414     {
    415       const Int_t   nspix = fSectorValid.At(sector);
    416       if (nspix == 0)
    417         continue;
    418      
    419       CalcSectorResults(fUsedEvents,nspix,sector);
    420     }
    421  
    422   fPedestalsOut->SetTotalEntries(fUsedEvents*fExtractWinSize);
    423   fPedestalsOut->SetReadyToSave();
    424 
    425   return;
     382Int_t MPedCalcPedRun::Finalize()
     383{
     384    if (fUsedEvents == 0)
     385        return kTRUE;
     386
     387    MRawEvtPixelIter pixel(fRawEvt);
     388    while (pixel.Next())
     389        CalcPixResults(fUsedEvents, pixel.GetPixelId());
     390
     391    //
     392    // Loop over the (two) area indices to get the averaged pedestal per aidx
     393    //
     394    for (UInt_t aidx=0; aidx<fAreaValid.GetSize(); aidx++)
     395        if (fAreaValid[aidx]>0)
     396            CalcAreaResults(fUsedEvents, fAreaValid[aidx], aidx);
     397
     398    //
     399    // Loop over the (six) sector indices to get the averaged pedestal per sector
     400    //
     401    for (UInt_t sector=0; sector<fSectorValid.GetSize(); sector++)
     402        if (fSectorValid[sector]>0)
     403            CalcSectorResults(fUsedEvents, fSectorValid[sector], sector);
     404
     405    fPedestalsOut->SetTotalEntries(fUsedEvents*fExtractWinSize);
     406    fPedestalsOut->SetReadyToSave();
     407
     408    return kTRUE;
     409}
     410
     411Int_t MPedCalcPedRun::PostProcess()
     412{
     413    if (!Finalize())
     414        return kFALSE;
     415
     416    return MExtractPedestal::PostProcess();
    426417}
    427418
     
    435426Bool_t MPedCalcPedRun::IsPedBitSet()
    436427{
    437 
    438   if (fRunHeader->GetRunNumber() < gkFirstRunWithFinalBits)
    439     return kFALSE;
    440      
    441   return (fEvtHeader->GetTriggerID() >> 3 & 1);
     428    if (fRunHeader->GetRunNumber() == 38996)
     429        return kTRUE;
     430
     431    if (fRunHeader->GetRunNumber() < gkFirstRunWithFinalBits)
     432        return kTRUE;
     433
     434    return (fEvtHeader->GetTriggerID() & BIT(3)) ? kTRUE : kFALSE;
    442435}
    443436
    444437void MPedCalcPedRun::Print(Option_t *o) const
    445438{
    446 
    447   MExtractPedestal::Print(o);
    448 
    449   *fLog << "Number overlap low-gain slices:           " << fOverlap << endl;
    450   *fLog << "First run out of sequence:                " << (fFirstRun?"yes":"no") << endl;
    451   *fLog << "Skip this run:                            " << (fSkip?"yes":"no") << endl;
    452   *fLog << "Number of used events so far:             " << fUsedEvents << endl;
    453   *fLog << endl;
    454 }
     439    *fLog << GetDescriptor() << ":" << endl;
     440    *fLog << "Name of input  MPedestalCam:   " << (fPedestalsIn?fPedestalsIn->GetName():fNamePedestalCamIn.Data()) << " (" << fPedestalsIn << ")" << endl;
     441    *fLog << "Name of output MPedestalCam:   " << (fPedestalsOut?fPedestalsOut->GetName():fNamePedestalCamOut.Data()) << " (" << fPedestalsOut << ")" << endl;
     442    *fLog << "ExtractWindow from slice       " << fExtractWinFirst << " to " << fExtractWinLast << " incl." << endl;
     443    *fLog << "Number overlap lo-gain slices: " << fOverlap << endl;
     444    *fLog << "First ped run out of sequence: " << (fIsFirstPedRun?"yes":"no") << endl;
     445    *fLog << "Number of used events so far:  " << fUsedEvents << endl;
     446
     447    if (fExtractor)
     448        *fLog << "Extractor used:                " << fExtractor->ClassName() << endl;
     449
     450    *fLog << endl;
     451}
  • trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.h

    r5550 r5558  
    66#endif
    77
    8 #ifndef MARS_MArrayD
    9 #include <MArrayD.h>
    10 #endif
    11 
    12 #ifndef MARS_MArrayI
    13 #include <MArrayI.h>
    14 #endif
    15 
    168class MRawEvtPixelIter;
    179class MPedestalPix;
     10
    1811class MPedCalcPedRun : public MExtractPedestal
    1912{
     13private:
     14    static const UShort_t fgExtractWinFirst;  // First FADC slice Hi-Gain (currently set to: 3)
     15    static const UShort_t fgExtractWinSize;   // Extraction Size Hi-Gain (currently set to: 14)
     16    static const UInt_t   gkFirstRunWithFinalBits; // First Run with pedestal trigger bit at place 3
    2017
    21   static const UShort_t fgExtractWinFirst;  // First FADC slice Hi-Gain (currently set to: 3)
    22   static const UShort_t fgExtractWinSize;   // Extraction Size Hi-Gain (currently set to: 14)
    23   static const UInt_t   gkFirstRunWithFinalBits; // First Run with pedestal trigger bit at place 3
    24  
    25   UShort_t fOverlap;         // Number of overlapping slices from High-Gain to Low-Gain
     18    UShort_t fOverlap;         // Number of overlapping slices from High-Gain to Low-Gain
    2619
    27   Bool_t  fFirstRun;         // Flag to tell if the first run out of many is used
    28   Bool_t  fSkip;             // Flag to tell if the Process has to be skipped
    29   ULong_t fUsedEvents;       // Number of used (not skipped) events
    30  
    31   Bool_t IsPedBitSet();
    32  
    33   Bool_t ReInit     (MParList *pList);
    34   Int_t  Process    ();
    35   Int_t  PostProcess();
    36  
    37   void ResetArrays();
    38   void CalcSums   ( MRawEvtPixelIter *pixel, Float_t &sum, UInt_t &ab0, UInt_t &ab1);
    39   void CalcExtractor( MRawEvtPixelIter *pixel, Float_t &sum, MPedestalPix &ped);
    40  
     20    Bool_t  fIsFirstPedRun;    //! Flag to tell if the first run out of many is used
     21    ULong_t fUsedEvents;       // Number of used (not skipped) events
     22
     23    Bool_t IsPedBitSet();
     24
     25    Bool_t ReInit(MParList *pList);
     26    Int_t  PreProcess(MParList *pList);
     27    Int_t  Process();
     28    Int_t  PostProcess();
     29
     30    void CheckExtractionWindow();
     31    void CalcSums(const MRawEvtPixelIter &pixel, Float_t &sum, UInt_t &ab0, UInt_t &ab1);
     32    void CalcExtractor(const MRawEvtPixelIter &pixel, Float_t &sum, MPedestalPix &ped);
     33
    4134public:
     35    MPedCalcPedRun(const char *name=NULL, const char *title=NULL);
    4236
    43   MPedCalcPedRun(const char *name=NULL, const char *title=NULL);
    44  
    45   void Finalize(); 
     37    void Print(Option_t *o="") const;
     38    Int_t Finalize();
    4639
    47   void Print(Option_t *o="") const;
    48  
    49   ClassDef(MPedCalcPedRun, 1)   // Task to calculate pedestals from pedestal runs
     40    ClassDef(MPedCalcPedRun, 1)   // Task to calculate pedestals from pedestal runs
    5041};
    5142
  • trunk/MagicSoft/Mars/mpedestal/MPedestalCam.cc

    r5516 r5558  
    108108  MPedestalCam &cam = (MPedestalCam&)obj;
    109109 
    110   const Int_t n = GetSize();
     110  Int_t n = GetSize();
    111111 
    112112  if (n==0)
     
    116116  for (int i=0; i<n; i++)
    117117    (*this)[i].Copy(cam[i]);
     118
     119  n = GetNumAverageArea();
     120  cam.InitAverageAreas(n);
     121  for (int i=0; i<n; i++)
     122    GetAverageArea(i).Copy(cam.GetAverageArea(i));
     123
     124  n = GetNumAverageSector();
     125  cam.InitAverageSectors(n);
     126  for (int i=0; i<n; i++)
     127    GetAverageSector(i).Copy(cam.GetAverageSector(i));
    118128}
    119129
     
    178188// independently if the MPedestalPix is filled with values or not.
    179189//
    180 const Int_t MPedestalCam::GetAverageAreas() const
     190const Int_t MPedestalCam::GetNumAverageArea() const
    181191{
    182192  return fAverageAreas->GetEntriesFast();
     
    188198// independently if the MPedestalPix is filled with values or not.
    189199//
    190 const Int_t MPedestalCam::GetAverageSectors() const
     200const Int_t MPedestalCam::GetNumAverageSector() const
    191201{
    192202  return fAverageSectors->GetEntriesFast();
     
    259269void MPedestalCam::Clear(Option_t *o)
    260270{
    261   fArray->ForEach(TObject, Clear)();
     271    { fArray->ForEach(TObject, Clear)(); }
     272    { fAverageAreas->ForEach(TObject, Clear)(); }
     273    { fAverageSectors->ForEach(TObject, Clear)(); }
    262274 
    263   //
    264   // another ForEach does not compile, thus have to do the loop ourselves:
    265   //
    266   for (Int_t i=0;i<GetAverageAreas();i++)
    267     fAverageAreas[i].Clear();
    268 
    269 
    270   //
    271   // another ForEach does not compile, thus have to do the loop ourselves:
    272   //
    273   for (Int_t i=0;i<GetAverageSectors();i++)
    274     fAverageSectors[i].Clear();
    275  
    276   fTotalEntries = 0;
     275    fTotalEntries = 0;
    277276}
    278277
     
    352351// arr[1]: Error (rms) of averaged pedestal (default: 0.)
    353352//
    354 // ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
    355 //
    356 TArrayF *MPedestalCam::GetAveragedPedPerArea(const MGeomCam &geom, const UInt_t ai, MBadPixelsCam *bad)
     353TArrayF MPedestalCam::GetAveragedPedPerArea(const MGeomCam &geom, const UInt_t ai, MBadPixelsCam *bad)
    357354{
    358355
     
    381378    }
    382379
    383   TArrayF *arr = new TArrayF(2);
    384   arr->AddAt(nr   ? mean/nr : -1.,0);
    385   arr->AddAt(nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0. ,1);
    386 
     380  TArrayF arr(2);
     381  arr[0] = nr   ? mean/nr : -1.;
     382  arr[1] = nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0;
    387383  return arr;
    388384}
     
    400396// arr[1]: Error (rms) of averaged pedestal RMS (default: 0.)
    401397//
    402 // ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
    403 //
    404 TArrayF *MPedestalCam::GetAveragedRmsPerArea(const MGeomCam &geom, const UInt_t ai, MBadPixelsCam *bad)
     398TArrayF MPedestalCam::GetAveragedRmsPerArea(const MGeomCam &geom, const UInt_t ai, MBadPixelsCam *bad)
    405399{
    406400
     
    428422    }
    429423
    430   TArrayF *arr = new TArrayF(2);
    431   arr->AddAt(nr   ? rms/nr : -1.,0);
    432   arr->AddAt(nr>1 ? TMath::Sqrt((rms2 - rms*rms/nr)/(nr-1)) : 0. ,1);
    433 
     424  TArrayF arr(2);
     425  arr[0] = nr   ? rms/nr : -1;
     426  arr[1] = nr>1 ? TMath::Sqrt((rms2 - rms*rms/nr)/(nr-1)) : 0;
    434427  return arr;
    435428}
     
    447440// arr[1]: Error (rms) of averaged pedestal (default: 0.)
    448441//
    449 // ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
    450 //
    451 TArrayF *MPedestalCam::GetAveragedPedPerSector(const MGeomCam &geom, const UInt_t sec, MBadPixelsCam *bad)
     442TArrayF MPedestalCam::GetAveragedPedPerSector(const MGeomCam &geom, const UInt_t sec, MBadPixelsCam *bad)
    452443{
    453444
     
    476467    }
    477468
    478   TArrayF *arr = new TArrayF(2);
    479   arr->AddAt(nr   ? mean/nr : -1.,0);
    480   arr->AddAt(nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0. ,1);
    481 
     469  TArrayF arr(2);
     470  arr[0] = nr   ? mean/nr : -1;
     471  arr[1] = nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0;
    482472  return arr;
    483473}
     
    495485// arr[1]: Error (rms) of averaged pedestal RMS (default: 0.)
    496486//
    497 // ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
    498 //
    499 TArrayF *MPedestalCam::GetAveragedRmsPerSector(const MGeomCam &geom, const UInt_t sec, MBadPixelsCam *bad)
     487TArrayF MPedestalCam::GetAveragedRmsPerSector(const MGeomCam &geom, const UInt_t sec, MBadPixelsCam *bad)
    500488{
    501489
     
    524512    }
    525513
    526   TArrayF *arr = new TArrayF(2);
    527   arr->AddAt(nr   ? rms/nr : -1.,0);
    528   arr->AddAt(nr>1 ? TMath::Sqrt((rms2 - rms*rms/nr)/(nr-1)) : 0. ,1);
    529 
     514  TArrayF arr(2);
     515  arr[0] = nr   ? rms/nr : -1;
     516  arr[1] = nr>1 ? TMath::Sqrt((rms2 - rms*rms/nr)/(nr-1)) : 0;
    530517  return arr;
    531518 
  • trunk/MagicSoft/Mars/mpedestal/MPedestalCam.h

    r5516 r5558  
    3636        MPedestalPix &GetAverageArea   ( UInt_t i );
    3737  const MPedestalPix &GetAverageArea   ( UInt_t i )            const;
    38   const Int_t         GetAverageAreas  ()                      const;
     38  const Int_t         GetNumAverageArea()                      const;
    3939        MPedestalPix &GetAverageSector ( UInt_t i );
    4040  const MPedestalPix &GetAverageSector ( UInt_t i )            const;
    41   const Int_t         GetAverageSectors()                      const;
     41  const Int_t         GetNumAverageSector()                    const;
    4242  Float_t             GetPedestalMin   ( const MGeomCam *cam ) const;
    4343  Float_t             GetPedestalMax   ( const MGeomCam *cam ) const;
     
    4545  ULong_t             GetTotalEntries  ()                      const { return fTotalEntries; }
    4646
    47   TArrayF *GetAveragedPedPerArea  ( const MGeomCam &geom, const UInt_t ai=0,  MBadPixelsCam *bad=NULL );
    48   TArrayF *GetAveragedPedPerSector( const MGeomCam &geom, const UInt_t sec=0, MBadPixelsCam *bad=NULL ); 
    49   TArrayF *GetAveragedRmsPerArea  ( const MGeomCam &geom, const UInt_t ai=0,  MBadPixelsCam *bad=NULL );
    50   TArrayF *GetAveragedRmsPerSector( const MGeomCam &geom, const UInt_t sec=0, MBadPixelsCam *bad=NULL ); 
     47  TArrayF GetAveragedPedPerArea  ( const MGeomCam &geom, const UInt_t ai=0,  MBadPixelsCam *bad=NULL );
     48  TArrayF GetAveragedPedPerSector( const MGeomCam &geom, const UInt_t sec=0, MBadPixelsCam *bad=NULL );
     49  TArrayF GetAveragedRmsPerArea  ( const MGeomCam &geom, const UInt_t ai=0,  MBadPixelsCam *bad=NULL );
     50  TArrayF GetAveragedRmsPerSector( const MGeomCam &geom, const UInt_t sec=0, MBadPixelsCam *bad=NULL );
    5151 
    5252        MPedestalPix &operator[]       ( Int_t i             );
  • trunk/MagicSoft/Mars/msignal/MExtractTime.cc

    r5548 r5558  
    215215void MExtractTime::Print(Option_t *o) const
    216216{
    217     *fLog << all;
    218217    if (IsA()==MExtractTime::Class())
    219218        *fLog << GetDescriptor() << ":" << endl;
  • trunk/MagicSoft/Mars/msignal/MExtractTimeAndCharge.cc

    r5548 r5558  
    115115  if (!fSignals)
    116116      return kFALSE;
     117
     118  *fLog << flush << inf;
     119  Print();
    117120 
    118121  return kTRUE;
     
    141144    fSignals->SetUsedFADCSlices(fHiGainFirst, fHiGainLast+fHiLoLast, fNumHiGainSamples,
    142145                                fLoGainFirst, fLoGainLast, fNumLoGainSamples);
    143  
    144   Print();
    145146
    146147  return kTRUE;
     
    244245void MExtractTimeAndCharge::Print(Option_t *o) const
    245246{
    246 
    247   *fLog << all;
    248247  if (IsA()==MExtractTimeAndCharge::Class())
    249248    *fLog << GetDescriptor() << ":" << endl;
    250249 
    251250  *fLog << dec << endl;
    252   *fLog << inf << "Taking " << fNumHiGainSamples
     251  *fLog << "Taking " << fNumHiGainSamples
    253252        << " HiGain samples from slice " << (Int_t)fHiGainFirst
    254253        << " to " << (Int_t)(fHiGainLast+fHiLoLast) << " incl" << endl;
  • trunk/MagicSoft/Mars/msignal/MExtractTimeAndChargeDigitalFilter.cc

    r5542 r5558  
    10121012void MExtractTimeAndChargeDigitalFilter::Print(Option_t *o) const
    10131013{
    1014 
    1015   *fLog << all;
    1016   *fLog << GetDescriptor() << ":" << endl;
    1017 
    1018   MExtractTimeAndCharge::Print(o);
    1019   *fLog << " Time Shift HiGain:  " << fTimeShiftHiGain << endl;
    1020   *fLog << " Time Shift LoGain:  " << fTimeShiftLoGain << endl;
    1021   *fLog << " Window Size HiGain: " << fWindowSizeHiGain << endl;
    1022   *fLog << " Window Size LoGain: " << fWindowSizeLoGain << endl;
    1023   *fLog << " Binning Res HiGain: " << fBinningResolutionHiGain << endl;
    1024   *fLog << " Binning Res LoGain: " << fBinningResolutionHiGain << endl;
    1025   *fLog << " Weights File:       " << fNameWeightsFile.Data() << endl;
    1026  
    1027   TString opt(o);
    1028   if (!opt.Contains("weights"))
    1029     return;
    1030  
    1031   *fLog << endl;
    1032   *fLog << inf << "Using the following weights: " << endl;
    1033   *fLog << "Hi-Gain:" << endl;
    1034   for (Int_t i=0; i<fBinningResolutionHiGain*fWindowSizeHiGain; i++)
    1035     *fLog << " " << fAmpWeightsHiGain[i] << " \t " << fTimeWeightsHiGain[i] << endl;
    1036  
    1037   *fLog << "Lo-Gain:" << endl;
    1038   for (Int_t i=0; i<fBinningResolutionLoGain*fWindowSizeLoGain; i++)
    1039     *fLog << " " << fAmpWeightsLoGain[i] << " \t " << fTimeWeightsLoGain[i] << endl;
    1040  
     1014    MExtractTimeAndCharge::Print(o);
     1015    *fLog << " Time Shift  HiGain: " << fTimeShiftHiGain         << "  LoGain: " << fTimeShiftLoGain << endl;
     1016    *fLog << " Window Size HiGain: " << fWindowSizeHiGain        << "  LoGain: " << fWindowSizeLoGain << endl;
     1017    *fLog << " Binning Res HiGain: " << fBinningResolutionHiGain << "  LoGain: " << fBinningResolutionHiGain << endl;
     1018    *fLog << " Weights File:       " << fNameWeightsFile.Data() << endl;
     1019
     1020    TString opt(o);
     1021    if (!opt.Contains("weights"))
     1022        return;
     1023
     1024    *fLog << endl;
     1025    *fLog << inf << "Using the following weights: " << endl;
     1026    *fLog << "Hi-Gain:" << endl;
     1027    for (Int_t i=0; i<fBinningResolutionHiGain*fWindowSizeHiGain; i++)
     1028        *fLog << " " << fAmpWeightsHiGain[i] << " \t " << fTimeWeightsHiGain[i] << endl;
     1029
     1030    *fLog << "Lo-Gain:" << endl;
     1031    for (Int_t i=0; i<fBinningResolutionLoGain*fWindowSizeLoGain; i++)
     1032        *fLog << " " << fAmpWeightsLoGain[i] << " \t " << fTimeWeightsLoGain[i] << endl;
    10411033}
  • trunk/MagicSoft/Mars/msignal/MExtractor.cc

    r5548 r5558  
    412412void MExtractor::Print(Option_t *o) const
    413413{
    414     *fLog << all;
    415 
    416414    if (IsA()==MExtractor::Class())
    417415        *fLog << GetDescriptor() << ":" << endl;
     
    419417    *fLog << " Hi Gain Range:  " << (int)fHiGainFirst << " " << (int)fHiGainLast << endl;
    420418    *fLog << " Lo Gain Range:  " << (int)fLoGainFirst << " " << (int)fLoGainLast << endl;
    421     *fLog << " Num Samples Hi: " << fNumHiGainSamples << endl;
    422     *fLog << " Num Samples Lo: " << fNumHiGainSamples << endl;
    423419    *fLog << " Saturation Lim: " << (int)fSaturationLimit << endl;
    424 }
     420    *fLog << " Num Samples HiGain: " << fNumHiGainSamples << "  LoGain: " << fNumLoGainSamples << endl;
     421}
Note: See TracChangeset for help on using the changeset viewer.