Changeset 4609 for trunk


Ignore:
Timestamp:
08/13/04 14:59:17 (20 years ago)
Author:
tbretz
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars
Files:
22 edited

Legend:

Unmodified
Added
Removed
  • trunk/MagicSoft/Mars/macros/readraw.C

    r948 r4609  
    1616!
    1717!
    18 !   Author(s): Thomas Bretz  12/2000 (tbretz@uni-sw.gwdg.de)
     18!   Author(s): Thomas Bretz  12/2000 <mailto:tbretz@astro.uni-wuerzburg.de>
    1919!
    20 !   Copyright: MAGIC Software Development, 2000-2001
     20!   Copyright: MAGIC Software Development, 2000-2004
    2121!
    2222!
     
    2424
    2525
    26 void readraw()
     26void readraw(const char *fname="/data/MAGIC/Period016/mcdata/spot_1cm/standard/gamma/Gamma_zbin9_90_7_1740to1749_w0.root")
    2727{
    2828    //
    2929    // open the file
    3030    //
    31     TFile input("test.root", "READ");
     31    TFile input(fname, "READ");
    3232
    3333    //
     
    4646    //
    4747    MRawRunHeader *runheader = new MRawRunHeader();
    48     runtree->GetBranch("MRawRunHeader")->SetAddress(&runheader);
     48    runtree->GetBranch("MRawRunHeader.")->SetAddress(&runheader);
    4949    runtree->GetEvent(0);
    5050    runheader->Print();
     
    5353    // open the Data Tree
    5454    //
    55     TTree *evttree = (TTree*) input.Get("Data") ;
     55    TTree *evttree = (TTree*) input.Get("Events");
    5656
    5757    //
    5858    // create the instances of the data to read in
    5959    //
    60     MRawEvtHeader  *evtheader = new MRawEvtHeader();
    61     MTime          *evttime   = new MTime();
    62     MRawEvtData    *evtdata   = new MRawEvtData();
    63     MRawCrateArray *evtcrate  = new MRawCrateArray();
     60    MRawEvtHeader  *evtheader = 0;
     61    MTime          *evttime   = 0;
     62    MRawEvtData    *evtdata   = 0;
     63    MRawCrateArray *evtcrate  = 0;
    6464
    6565    //
    6666    // enable the corresponding branches
    6767    //
    68     evttree->GetBranch("MRawEvtHeader")->SetAddress(&evtheader);
    69     evttree->GetBranch("MTime")->SetAddress(&evttime);
    70     evttree->GetBranch("MRawEvtData")->SetAddress(&evtdata);
    71     evttree->GetBranch("MRawCrateArray")->SetAddress(&evtcrate);
     68    evttree->GetBranch("MRawEvtHeader.")->SetAddress(&evtheader);
     69    evttree->GetBranch("MRawEvtData.")->SetAddress(&evtdata);
     70
     71    // Use this for real data only
     72    //evttree->GetBranch("MTime.")->SetAddress(&evttime);
     73    //evttree->GetBranch("MRawCrateArray.")->SetAddress(&evtcrate);
    7274
    7375    //
     
    8385
    8486      evtheader->Print();
    85       evttime->Print();
    86       evtcrate->Print();
    8787      evtdata->Print();
     88
     89      // Use this for real data only
     90      //evttime->Print();
     91      //evtcrate->Print();
    8892    }
    8993}
  • trunk/MagicSoft/Mars/mbase/MTaskEnv.cc

    r4601 r4609  
    172172    if (!fTask)
    173173    {
    174         *fLog << err << GetDescriptor() << " - ERROR: Not task setup." << endl;
     174        *fLog << err << GetDescriptor() << " - ERROR: No task setup." << endl;
    175175        return kFALSE;
    176176    }
  • trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.cc

    r4603 r4609  
    15241524        }
    15251525
    1526       *fLog << inf << GetDescriptor() << ": Mean overall F-Factor (photons to FADC counts) "
    1527             << "with area index: " << i << ": "
    1528             << Form("%4.2f+-%4.2f",mean,sigma) << endl;
     1526      *fLog << inf << endl;
     1527      *fLog << inf << GetDescriptor() << ": Mean F-Factor "
     1528          << "with area index #" << i << ": "
     1529            << Form("%4.2f+-%4.2f",mean,sigma) << "fadc/ph" << endl;
    15291530
    15301531      lowlim  [i] = 1.1;
     
    15501551      if ( ffactor < lowlim[aidx] || ffactor > upplim[aidx] )
    15511552        {
    1552           *fLog << warn << GetDescriptor() << ": Overall F-Factor: "
    1553                 << Form("%5.2f out of %3.1f sigma limit: ",ffactor,fFFactorErrLimit)
    1554                 << Form("[%5.2f,%5.2f] pixel %4i",lowlim[aidx],upplim[aidx],i) << endl;
     1553          *fLog << warn << GetDescriptor() << ": Overall F-Factor "
     1554                << Form("%5.2f",ffactor) << " out of range ["
     1555                << Form("%5.2f,%5.2f",lowlim[aidx],upplim[aidx]) << "] pixel " << i << endl;
     1556
    15551557          bad.SetUncalibrated( MBadPixelsPix::kDeviatingFFactor );
    15561558          bad.SetUnsuitable  ( MBadPixelsPix::kUnsuitableRun    );
     
    18571859  return Form("%s/%s", (const char*)fOutputPath, (const char*)fOutputFile);
    18581860}
     1861
     1862Int_t MCalibrationChargeCalc::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
     1863{
     1864    Bool_t rc = kFALSE;
     1865    if (IsEnvDefined(env, prefix, "ChargeLimit", print))
     1866    {
     1867        SetChargeLimit(GetEnvValue(env, prefix, "ChargeLimit", fChargeLimit));
     1868        rc = kTRUE;
     1869    }
     1870
     1871    if (IsEnvDefined(env, prefix, "ChargeErrLimit", print))
     1872    {
     1873        SetChargeErrLimit(GetEnvValue(env, prefix, "ChargeErrLimit", fChargeErrLimit));
     1874        rc = kTRUE;
     1875    }
     1876
     1877    if (IsEnvDefined(env, prefix, "ChargeRelErrLimit", print))
     1878    {
     1879        SetChargeRelErrLimit(GetEnvValue(env, prefix, "ChargeRelErrLimit", fChargeRelErrLimit));
     1880        rc = kTRUE;
     1881    }
     1882    if (IsEnvDefined(env, prefix, "Debug", print))
     1883    {
     1884        SetDebug(GetEnvValue(env, prefix, "Debug", IsDebug()));
     1885        rc = kTRUE;
     1886    }
     1887    if (IsEnvDefined(env, prefix, "FFactorErrLimit", print))
     1888    {
     1889        SetFFactorErrLimit(GetEnvValue(env, prefix, "FFactorErrLimit", fFFactorErrLimit));
     1890        rc = kTRUE;
     1891    }
     1892    if (IsEnvDefined(env, prefix, "LambdaErrLimit", print))
     1893    {
     1894        SetLambdaErrLimit(GetEnvValue(env, prefix, "LambdaErrLimit", fLambdaErrLimit));
     1895        rc = kTRUE;
     1896    }
     1897    if (IsEnvDefined(env, prefix, "LambdaCheckLimit", print))
     1898    {
     1899        SetLambdaCheckLimit(GetEnvValue(env, prefix, "LambdaCheckLimit", fLambdaCheckLimit));
     1900        rc = kTRUE;
     1901    }
     1902    if (IsEnvDefined(env, prefix, "PheErrLimit", print))
     1903    {
     1904        SetPheErrLimit(GetEnvValue(env, prefix, "PheErrLimit", fPheErrLimit));
     1905        rc = kTRUE;
     1906    }
     1907    // void SetPulserColor(const MCalibrationCam::PulserColor_t col) { fPulserColor = col; }
     1908
     1909    return rc;
     1910}
  • trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.h

    r4603 r4609  
    127127  Int_t  Process    ();
    128128  Int_t  PostProcess();
     129  Int_t  ReadEnv(const TEnv &env, TString prefix, Bool_t print);
    129130
    130131public:
  • trunk/MagicSoft/Mars/mcalib/MHCalibrationCam.cc

    r4603 r4609  
    10231023}
    10241024
     1025Int_t MHCalibrationCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
     1026{
     1027    Bool_t rc = kFALSE;
     1028    if (IsEnvDefined(env, prefix, "Debug", print))
     1029    {
     1030        SetDebug(GetEnvValue(env, prefix, "Debug", IsDebug()));
     1031        rc = kTRUE;
     1032    }
     1033
     1034    return rc;
     1035}
  • trunk/MagicSoft/Mars/mcalib/MHCalibrationCam.h

    r4602 r4609  
    109109  void InitHists(MHGausEvents &hist, MBadPixelsPix &bad, const Int_t i);
    110110
     111  Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print);
     112
    111113public:
    112114
  • trunk/MagicSoft/Mars/mjobs/MJCalibration.cc

    r4607 r4609  
    8787#include "MJCalibration.h"
    8888
     89#include <TEnv.h>
    8990#include <TFile.h>
    9091#include <TStyle.h>
     
    9697
    9798#include "MRunIter.h"
     99#include "MSequence.h"
    98100#include "MParList.h"
    99101#include "MTaskList.h"
     
    125127#include "MRawFileRead.h"
    126128#include "MGeomApply.h"
     129#include "MTaskEnv.h"
    127130#include "MBadPixelsMerge.h"
    128131#include "MBadPixelsCam.h"
     
    158161//
    159162MJCalibration::MJCalibration(const char *name, const char *title)
    160     : fRuns(0), fExtractor(NULL), fTimeExtractor(NULL),
     163    : fEnv(0), fRuns(0), fSequence(0), fExtractor(NULL), fTimeExtractor(NULL),
    161164      fColor(MCalibrationCam::kNONE), fDisplayType(kNormalDisplay),
    162165      fRelTimes(kFALSE), fDataCheck(kFALSE), fDebug(kFALSE)
     
    171174}
    172175
     176MJCalibration::~MJCalibration()
     177{
     178    if (fEnv)
     179        delete fEnv;
     180}
    173181
    174182// --------------------------------------------------------------------------
     
    224232    TString title = fDisplay->GetTitle();
    225233    title += "--  Calibration ";
    226     title += fRuns->GetRunsAsString();
     234    title += fSequence ? Form("calib%06d", fSequence->GetSequence()) : fRuns->GetRunsAsString();
    227235    title += "  --";
    228236    fDisplay->SetTitle(title);
     
    234242
    235243    // Create histograms to display
    236     MHCamera disp1 (geomcam, Form("%s%s","Charge",(fRuns->GetRunsAsFileName()).Data()),
    237                     "Fitted Mean Charges");
    238     MHCamera disp2 (geomcam, Form("%s%s","SigmaCharge",(fRuns->GetRunsAsFileName()).Data()),
    239                     "Sigma of Fitted Charges");
    240     MHCamera disp3 (geomcam, Form("%s%s","RSigma",(fRuns->GetRunsAsFileName()).Data()),
    241                     "Reduced Sigmas");
    242     MHCamera disp4 (geomcam, Form("%s%s","RSigmaPerCharge",(fRuns->GetRunsAsFileName()).Data()), 
    243                     "Reduced Sigma per Charge");
    244     MHCamera disp5 (geomcam, Form("%s%s","NumPhes",(fRuns->GetRunsAsFileName()).Data()),
    245                     "Nr. of Phe's (F-Factor Method)");
    246     MHCamera disp6 (geomcam, Form("%s%s","ConvFADC2Phes",(fRuns->GetRunsAsFileName()).Data()),
    247                     "Conversion Factor (F-Factor Method)");
    248     MHCamera disp7 (geomcam, Form("%s%s","TotalFFactor",(fRuns->GetRunsAsFileName()).Data()),
    249                     "Total F-Factor (F-Factor Method)");
    250     MHCamera disp8 (geomcam, Form("%s%s","CascadesQEFFactor",(fRuns->GetRunsAsFileName()).Data()),
    251                     "Cascades QE (F-Factor Method)");
    252     MHCamera disp9 (geomcam, Form("%s%s","CascadesQEBlindPix",(fRuns->GetRunsAsFileName()).Data()),
    253                     "Cascades QE (Blind Pixel Method)");
    254     MHCamera disp10(geomcam, Form("%s%s","CascadesQEPINDiode",(fRuns->GetRunsAsFileName()).Data()),
    255                     "Cascades QE (PIN Diode Method)");
    256     MHCamera disp11(geomcam, Form("%s%s","CascadesQECombined",(fRuns->GetRunsAsFileName()).Data()),
    257                     "Cascades QE (Combined Method)");
    258     MHCamera disp12(geomcam, Form("%s%s","FFactorValid",(fRuns->GetRunsAsFileName()).Data()),
    259                     "Pixels with valid F-Factor calibration");
    260     MHCamera disp13(geomcam, Form("%s%s","BlindPixelValid",(fRuns->GetRunsAsFileName()).Data()),
    261                     "Pixels with valid BlindPixel calibration");
    262     MHCamera disp14(geomcam, Form("%s%s","PINdiodeValid",(fRuns->GetRunsAsFileName()).Data()),
    263                     "Pixels with valid PINDiode calibration");
    264     MHCamera disp15(geomcam, Form("%s%s","CombinedValid",(fRuns->GetRunsAsFileName()).Data()),
    265                     "Pixels with valid Combined calibration");
    266     MHCamera disp16(geomcam, Form("%s%s","Saturation",(fRuns->GetRunsAsFileName()).Data()),
    267                     "Pixels with saturated Hi Gain");
    268     MHCamera disp17(geomcam, Form("%s%s","ConversionMeans",(fRuns->GetRunsAsFileName()).Data()),
    269                     "Conversion HiGain.vs.LoGain Means");
    270     MHCamera disp18(geomcam, Form("%s%s","ConversionSigmas",(fRuns->GetRunsAsFileName()).Data()),
    271                     "Conversion HiGain.vs.LoGain Sigmas");
    272     MHCamera disp19(geomcam, Form("%s%s","HiGainPickup",(fRuns->GetRunsAsFileName()).Data()),
    273                     "Number Pickup events Hi Gain");
    274     MHCamera disp20(geomcam, Form("%s%s","LoGainPickup",(fRuns->GetRunsAsFileName()).Data()),
    275                     "Number Pickup events Lo Gain");
    276     MHCamera disp21(geomcam, Form("%s%s","HiGainBlackout",(fRuns->GetRunsAsFileName()).Data()),
    277                     "Number Blackout events Hi Gain");
    278     MHCamera disp22(geomcam, Form("%s%s","LoGainBlackout",(fRuns->GetRunsAsFileName()).Data()),
    279                     "Number Blackout events Lo Gain");
    280     MHCamera disp23(geomcam, Form("%s%s","Excluded",(fRuns->GetRunsAsFileName()).Data()),
    281                     "Pixels previously excluded");
    282     MHCamera disp24(geomcam, Form("%s%s","UnSuitable",(fRuns->GetRunsAsFileName()).Data()),
    283                     "Pixels not suited for further analysis");
    284     MHCamera disp25(geomcam, Form("%s%s","UnReliable",(fRuns->GetRunsAsFileName()).Data()),
    285                     "Pixels not reliable for further analysis");
    286     MHCamera disp26(geomcam, Form("%s%s","HiGainOscillating",(fRuns->GetRunsAsFileName()).Data()),
    287                     "Oscillating Pixels High Gain");
    288     MHCamera disp27(geomcam, Form("%s%s","LoGainOscillating",(fRuns->GetRunsAsFileName()).Data()),
    289                     "Oscillating Pixels Low Gain");
    290     MHCamera disp28(geomcam, Form("%s%s","AbsTimeMean",(fRuns->GetRunsAsFileName()).Data()),
    291                     "Abs. Arrival Times");
    292     MHCamera disp29(geomcam, Form("%s%s","AbsTimeRms",(fRuns->GetRunsAsFileName()).Data()),
    293                     "RMS of Arrival Times");
    294     MHCamera disp30(geomcam, Form("%s%s","MeanTime",(fRuns->GetRunsAsFileName()).Data()),
    295                     "Mean Rel. Arrival Times");
    296     MHCamera disp31(geomcam, Form("%s%s","SigmaTime",(fRuns->GetRunsAsFileName()).Data()),
    297                     "Sigma Rel. Arrival Times");
    298     MHCamera disp32(geomcam, Form("%s%s","TimeProb",(fRuns->GetRunsAsFileName()).Data()),
    299                     "Probability of Time Fit");
    300     MHCamera disp33(geomcam, Form("%s%s","TimeNotFitValid",(fRuns->GetRunsAsFileName()).Data()),
    301                     "Pixels with not valid fit results");
    302     MHCamera disp34(geomcam, Form("%s%s","TimeOscillating",(fRuns->GetRunsAsFileName()).Data()),
    303                     "Oscillating Pixels");
     244    MHCamera disp1 (geomcam, "Charge",            "Fitted Mean Charges");
     245    MHCamera disp2 (geomcam, "SigmaCharge",       "Sigma of Fitted Charges");
     246    MHCamera disp3 (geomcam, "RSigma",            "Reduced Sigmas");
     247    MHCamera disp4 (geomcam, "RSigmaPerCharge",   "Reduced Sigma per Charge");
     248    MHCamera disp5 (geomcam, "NumPhes",           "Nr. of Phe's (F-Factor Method)");
     249    MHCamera disp6 (geomcam, "ConvFADC2Phes",     "Conversion Factor (F-Factor Method)");
     250    MHCamera disp7 (geomcam, "TotalFFactor",      "Total F-Factor (F-Factor Method)");
     251    MHCamera disp8 (geomcam, "CascadesQEFFactor", "Cascades QE (F-Factor Method)");
     252    MHCamera disp9 (geomcam, "CascadesQEBlindPix","Cascades QE (Blind Pixel Method)");
     253    MHCamera disp10(geomcam, "CascadesQEPINDiode","Cascades QE (PIN Diode Method)");
     254    MHCamera disp11(geomcam, "CascadesQECombined","Cascades QE (Combined Method)");
     255    MHCamera disp12(geomcam, "FFactorValid",      "Pixels with valid F-Factor calibration");
     256    MHCamera disp13(geomcam, "BlindPixelValid",   "Pixels with valid BlindPixel calibration");
     257    MHCamera disp14(geomcam, "PINdiodeValid",     "Pixels with valid PINDiode calibration");
     258    MHCamera disp15(geomcam, "CombinedValid",     "Pixels with valid Combined calibration");
     259    MHCamera disp16(geomcam, "Saturation",        "Pixels with saturated Hi Gain");
     260    MHCamera disp17(geomcam, "ConversionMeans",   "Conversion HiGain.vs.LoGain Means");
     261    MHCamera disp18(geomcam, "ConversionSigmas",  "Conversion HiGain.vs.LoGain Sigmas");
     262    MHCamera disp19(geomcam, "HiGainPickup",      "Number Pickup events Hi Gain");
     263    MHCamera disp20(geomcam, "LoGainPickup",      "Number Pickup events Lo Gain");
     264    MHCamera disp21(geomcam, "HiGainBlackout",    "Number Blackout events Hi Gain");
     265    MHCamera disp22(geomcam, "LoGainBlackout",    "Number Blackout events Lo Gain");
     266    MHCamera disp23(geomcam, "Excluded",          "Pixels previously excluded");
     267    MHCamera disp24(geomcam, "UnSuitable",        "Pixels not suited for further analysis");
     268    MHCamera disp25(geomcam, "UnReliable",        "Pixels not reliable for further analysis");
     269    MHCamera disp26(geomcam, "HiGainOscillating", "Oscillating Pixels High Gain");
     270    MHCamera disp27(geomcam, "LoGainOscillating", "Oscillating Pixels Low Gain");
     271    MHCamera disp28(geomcam, "AbsTimeMean",       "Abs. Arrival Times");
     272    MHCamera disp29(geomcam, "AbsTimeRms",        "RMS of Arrival Times");
     273    MHCamera disp30(geomcam, "MeanTime",          "Mean Rel. Arrival Times");
     274    MHCamera disp31(geomcam, "SigmaTime",         "Sigma Rel. Arrival Times");
     275    MHCamera disp32(geomcam, "TimeProb",          "Probability of Time Fit");
     276    MHCamera disp33(geomcam, "TimeNotFitValid",   "Pixels with not valid fit results");
     277    MHCamera disp34(geomcam, "TimeOscillating",   "Oscillating Pixels");
    304278
    305279    // Fitted charge means and sigmas
     
    665639Bool_t MJCalibration::FindColor()
    666640{
    667 
    668   const UInt_t nruns = fRuns->GetNumRuns();
     641    if (fSequence)
     642    {
     643        fColor = MCalibrationCam::kCT1;
     644        return kTRUE;
     645    }
     646
     647const UInt_t nruns = fRuns->GetNumRuns();
    669648
    670649  if (nruns == 0)
     
    880859    }
    881860
    882   return kTRUE;
    883 }
    884 
    885 
    886 
     861
     862    return kTRUE;
     863}
    887864
    888865// --------------------------------------------------------------------------
     
    892869const char* MJCalibration::GetOutputFile() const
    893870{
    894 
    895   if (!fRuns)
    896     return "";
    897  
    898   return Form("%s/%s-F1.root", (const char*)fOutputPath, (const char*)fRuns->GetRunsAsFileName());
    899 }
    900 
     871    if (fSequence)
     872        return Form("%s/calib%06d.root", (const char*)fOutputPath, fSequence->GetSequence());
     873    if (!fRuns)
     874        return "";
     875
     876    return Form("%s/%s-F1.root", (const char*)fOutputPath, (const char*)fRuns->GetRunsAsFileName());
     877}
    901878
    902879Bool_t MJCalibration::IsUseBlindPixel() const
    903880{
    904   return TESTBIT(fDevices,kUseBlindPixel);
    905 }
    906 
     881    return TESTBIT(fDevices,kUseBlindPixel);
     882}
    907883
    908884Bool_t MJCalibration::IsUsePINDiode() const
    909885{
    910   return TESTBIT(fDevices,kUsePINDiode);
    911 }
    912 
    913  
    914 
    915 
     886    return TESTBIT(fDevices,kUsePINDiode);
     887}
     888
     889void MJCalibration::SetEnv(const char *env)
     890{
     891    if (fEnv)
     892        delete fEnv;
     893    fEnv = new TEnv(env);
     894}
     895
     896void MJCalibration::CheckEnv()
     897{
     898    if (!fEnv)
     899        return;
     900
     901    TString e1 = fEnv->GetValue("MJCalibration.OutputPath", "");
     902    if (!e1.IsNull())
     903    {
     904        e1.ReplaceAll("\015", "");
     905        SetOutputPath(e1);
     906    }
     907}
    916908
    917909// --------------------------------------------------------------------------
     
    925917
    926918    return kTRUE;
     919}
     920
     921void MJCalibration::InitBlindPixel(MExtractBlindPixel &blindext,
     922                                   MHCalibrationChargeBlindCam &blindcam)
     923{
     924    Int_t run = fSequence ? fSequence->GetLastRun() : fRuns->GetRuns()[fRuns->GetNumRuns()-1];
     925
     926    //
     927    // Initialize the blind pixel. Unfortunately, there is a hardware difference
     928    // in the first blind pixel until run "gkSecondBlindPixelInstallation" and the
     929    // later setup. The first needs to use a filter because of the length of
     930    // spurious NSB photon signals. The latter get better along extracting the amplitude
     931    // from a small window.
     932    //
     933    if (run < gkSecondBlindPixelInstallation)
     934    {
     935        blindext.SetModified(kFALSE);
     936        blindext.SetExtractionType(MExtractBlindPixel::kIntegral);
     937        blindext.SetExtractionType(MExtractBlindPixel::kFilter);
     938        blindext.SetRange(10,19,0,6);
     939        blindext.SetNSBFilterLimit(70);
     940        blindcam.SetFitFunc( MHCalibrationChargeBlindPix::kEPoisson5 );
     941    }
     942    else
     943    {
     944        blindext.SetModified(kTRUE);
     945        blindext.SetExtractionType(MExtractBlindPixel::kAmplitude);
     946        blindext.SetExtractionType(MExtractBlindPixel::kFilter);
     947        blindext.SetRange(5,8,0,2);
     948        blindext.SetNSBFilterLimit(38);
     949
     950        if (run < gkThirdBlindPixelInstallation)
     951            blindext.SetNumBlindPixels(2);
     952        else
     953            blindext.SetNumBlindPixels(3);
     954    }
    927955}
    928956
     
    965993Bool_t MJCalibration::ProcessFile(MPedestalCam &pedcam)
    966994{
    967   if (!fRuns)
    968     {
    969       *fLog << err << "No Runs choosen... abort." << endl;
    970       return kFALSE;
    971     }
     995    if (!fRuns && !fSequence)
     996    {
     997        *fLog << err << "No Runs choosen... abort." << endl;
     998        return kFALSE;
     999    }
     1000
     1001    if (!fSequence && fRuns->GetNumRuns() != fRuns->GetNumEntries())
     1002    {
     1003        *fLog << err << "Number of files found doesn't match number of runs... abort."
     1004            << fRuns->GetNumRuns() << " vs. " << fRuns->GetNumEntries() << endl;
     1005        return kFALSE;
     1006    }
     1007
     1008    *fLog << inf;
     1009    fLog->Separator(GetDescriptor());
     1010
     1011    if (!FindColor())
     1012        return kFALSE;
     1013
     1014    *fLog << "Calculate MCalibrationCam from ";
     1015    if (fSequence)
     1016        *fLog << "Sequence #" << fSequence->GetSequence() << endl;
     1017    else
     1018        *fLog << "Runs " << fRuns->GetRunsAsString() << endl;
     1019    *fLog << endl;
     1020
     1021    // Setup Tasklist
     1022    MParList plist;
     1023    MTaskList tlist;
     1024    plist.AddToList(&tlist);
     1025
     1026    MReadMarsFile read("Events");
     1027    MRawFileRead rawread(NULL);
     1028
     1029    MDirIter iter;
     1030    if (fSequence)
     1031        fSequence->SetupCalRuns(iter);
     1032
     1033    if (fDataCheck)
     1034    {
     1035        rawread.AddFiles(fSequence ? iter : *fRuns);
     1036        tlist.AddToList(&rawread);
     1037    }
     1038    else
     1039    {
     1040        read.DisableAutoScheme();
     1041        static_cast<MRead&>(read).AddFiles(fSequence ? iter : *fRuns);
     1042        tlist.AddToList(&read);
     1043    }
     1044
     1045    MHCalibrationChargeCam      chargecam;
     1046    MHCalibrationChargeBlindCam blindcam;
     1047
     1048    plist.AddToList(&pedcam);
     1049    plist.AddToList(&chargecam);
     1050    plist.AddToList(&blindcam);
     1051    plist.AddToList(&fBadPixels);
     1052    plist.AddToList(&fQECam);
     1053    plist.AddToList(&fCalibrationCam);
     1054    plist.AddToList(&fCalibrationBlindCam);
     1055    plist.AddToList(&fCalibrationPINDiode);
     1056    plist.AddToList(&fRelTimeCam);
     1057
     1058    MGeomApply               apply;
     1059    MBadPixelsMerge          merge(&fBadPixels);
     1060    MExtractPINDiode         pinext;
     1061    MExtractBlindPixel       blindext;
     1062    InitBlindPixel(blindext, blindcam);
     1063    MExtractSlidingWindow    extract2;
     1064    MExtractTimeFastSpline   timespline;
     1065    MCalibrationChargeCalc   calcalc;
     1066    if (!fSequence)
     1067    {
     1068        calcalc.SetOutputPath(fOutputPath);
     1069        calcalc.SetOutputFile(Form("%s-ChargeCalibStat.txt",(const char*)fRuns->GetRunsAsFileName()));
     1070    }
     1071
     1072    if (fDebug)
     1073    {
     1074        chargecam.SetDebug();
     1075        calcalc.SetDebug();
     1076    }
     1077
     1078    //
     1079    // As long as there are no DM's, have to colour by hand
     1080    //
     1081    calcalc.SetPulserColor(fColor);
     1082
     1083    MFillH fillpin("MHCalibrationChargePINDiode", "MExtractedSignalPINDiode");
     1084    MFillH fillbnd("MHCalibrationChargeBlindPix", "MExtractedSignalBlindPixel");
     1085    MFillH fillcam("MHCalibrationChargeCam",      "MExtractedSignalCam");
     1086    MFillH filltme("MHCalibrationRelTimeCam",     "MArrivalTimeCam");
     1087    fillpin.SetNameTab("PINDiode");
     1088    fillbnd.SetNameTab("BlindPix");
     1089    fillcam.SetNameTab("Charge");
     1090    filltme.SetNameTab("RelTimes");
    9721091 
    973   if (fRuns->GetNumRuns() != fRuns->GetNumEntries())
    974     {
    975       *fLog << err << "Number of files found doesn't match number of runs... abort."
    976             << fRuns->GetNumRuns() << " vs. " << fRuns->GetNumEntries() << endl;
    977       return kFALSE;
    978     }
    979 
    980   *fLog << inf;
    981   fLog->Separator(GetDescriptor());
    982 
    983   if (!FindColor())
    984     return kFALSE;
    985  
    986   *fLog << "Calculate MCalibrationCam from Runs " << fRuns->GetRunsAsString() << endl;
    987   *fLog << endl;
    988  
    989   // Setup Tasklist
    990   MParList plist;
    991   MTaskList tlist;
    992   plist.AddToList(&tlist);
    993  
    994   MReadMarsFile read("Events");
    995   MRawFileRead rawread(NULL);
    996 
    997   if (fDataCheck)
    998   {
    999      rawread.AddFiles(*fRuns);
    1000      tlist.AddToList(&rawread);
    1001   }
    1002   else
    1003   {
    1004       read.DisableAutoScheme();
    1005       static_cast<MRead&>(read).AddFiles(*fRuns);
    1006       tlist.AddToList(&read);
    1007   }
    1008 
    1009   MHCalibrationChargeCam      chargecam;
    1010   MHCalibrationChargeBlindCam blindcam;
    1011 
    1012   plist.AddToList(&pedcam);
    1013   plist.AddToList(&chargecam); 
    1014   plist.AddToList(&blindcam); 
    1015   plist.AddToList(&fBadPixels);
    1016   plist.AddToList(&fQECam);
    1017   plist.AddToList(&fCalibrationCam);
    1018   plist.AddToList(&fCalibrationBlindCam);
    1019   plist.AddToList(&fCalibrationPINDiode);
    1020   plist.AddToList(&fRelTimeCam);
    1021 
    1022   MGeomApply               apply;
    1023   MExtractPINDiode         pinext;
    1024   MExtractBlindPixel       blindext;
    1025  
    1026   //
    1027   // Initialize the blind pixel. Unfortunately, there is a hardware difference
    1028   // in the first blind pixel until run "gkSecondBlindPixelInstallation" and the
    1029   // later setup. The first needs to use a filter because of the length of
    1030   // spurious NSB photon signals. The latter get better along extracting the amplitude
    1031   // from a small window.
    1032   //
    1033   TArrayI arr = fRuns->GetRuns();
    1034   if (arr[fRuns->GetNumRuns()-1] < gkSecondBlindPixelInstallation)
    1035     {
    1036       blindext.SetModified(kFALSE);
    1037       blindext.SetExtractionType(MExtractBlindPixel::kIntegral);
    1038       blindext.SetExtractionType(MExtractBlindPixel::kFilter);
    1039       blindext.SetRange(10,19,0,6);
    1040       blindext.SetNSBFilterLimit(70);
    1041       blindcam.SetFitFunc( MHCalibrationChargeBlindPix::kEPoisson5 );
    1042     }
    1043   else
    1044     {
    1045       blindext.SetModified(kTRUE);
    1046       blindext.SetExtractionType(MExtractBlindPixel::kAmplitude);
    1047       blindext.SetExtractionType(MExtractBlindPixel::kFilter);
    1048       blindext.SetRange(5,8,0,2);
    1049       blindext.SetNSBFilterLimit(38);
    1050       if (arr[fRuns->GetNumRuns()-1] < gkThirdBlindPixelInstallation)
    1051         blindext.SetNumBlindPixels(2);
    1052       else
    1053         blindext.SetNumBlindPixels(3);
    1054     }
    1055  
    1056   MExtractSlidingWindow    extract2;
    1057   MExtractTimeFastSpline   timespline;
    1058   MCalibrationChargeCalc   calcalc;
    1059   calcalc.SetOutputPath(fOutputPath);
    1060   calcalc.SetOutputFile(Form("%s-ChargeCalibStat.txt",(const char*)fRuns->GetRunsAsFileName()));
    1061 
    1062   if (fDebug)
    1063     {
    1064       chargecam.SetDebug();
    1065       calcalc.SetDebug();     
    1066     }
    1067 
    1068   MCalibrationRelTimeCalc  timecalc;
    1069   timecalc.SetOutputPath(fOutputPath);
    1070   timecalc.SetOutputFile(Form("%s-TimeCalibStat.txt",(const char*)fRuns->GetRunsAsFileName()));
    1071  
    1072   //
    1073   // As long as there are no DM's, have to colour by hand
    1074   //
    1075   calcalc.SetPulserColor(fColor);
    1076  
    1077   MFillH fillpin("MHCalibrationChargePINDiode", "MExtractedSignalPINDiode");
    1078   MFillH fillbnd("MHCalibrationChargeBlindCam", "MExtractedSignalBlindPixel");
    1079   MFillH fillcam("MHCalibrationChargeCam",      "MExtractedSignalCam");
    1080   MFillH filltme("MHCalibrationRelTimeCam",     "MArrivalTimeCam");
    1081   fillpin.SetNameTab("PINDiode");
    1082   fillbnd.SetNameTab("BlindPix");
    1083   fillcam.SetNameTab("Charge");
    1084   filltme.SetNameTab("RelTimes");
    1085 
    1086   TString drawoption;
    1087 
    1088   if (fDisplayType == kDataCheckDisplay)
    1089     drawoption += "datacheck";
    1090   if (fDisplayType == kFullDisplay)
    1091     drawoption += " all";
    1092  
    1093   fillcam.SetDrawOption(drawoption.Data());
    1094   fillbnd.SetDrawOption(drawoption.Data());
    1095   fillpin.SetDrawOption(drawoption.Data());
    1096   filltme.SetDrawOption(drawoption.Data());
    1097   //
    1098   // Apply a filter against cosmics
    1099   // (will have to be needed in the future
    1100   // when the calibration hardware-trigger is working)
    1101   //
    1102   MFCosmics cosmics;
    1103   MContinue cont(&cosmics);
    1104  
    1105   //    tlist.AddToList(&merge);
    1106   tlist.AddToList(&apply);
    1107 
    1108   if (fExtractor)
    1109     tlist.AddToList(fExtractor);
    1110   else
    1111     {
    1112       *fLog << warn << GetDescriptor()
    1113             << ": No extractor has been chosen, take default MExtractSlidingWindow " << endl;
    1114       tlist.AddToList(&extract2);
    1115     }
    1116  
    1117   tlist.AddToList(&pinext); 
    1118   tlist.AddToList(&blindext);
    1119  
    1120   if (fRelTimes)
    1121     {
    1122       if (fTimeExtractor)
    1123         tlist.AddToList(fTimeExtractor);
    1124       else
     1092    TString drawoption;
     1093
     1094    if (fDisplayType == kDataCheckDisplay)
     1095        drawoption += "datacheck";
     1096    if (fDisplayType == kFullDisplay)
     1097        drawoption += " all";
     1098
     1099    fillcam.SetDrawOption(drawoption.Data());
     1100    fillbnd.SetDrawOption(drawoption.Data());
     1101    fillpin.SetDrawOption(drawoption.Data());
     1102    filltme.SetDrawOption(drawoption.Data());
     1103
     1104    //
     1105    // Apply a filter against cosmics
     1106    // (will have to be needed in the future
     1107    // when the calibration hardware-trigger is working)
     1108    //
     1109    MFCosmics cosmics;
     1110    MContinue cont(&cosmics);
     1111
     1112    tlist.AddToList(&merge);
     1113    tlist.AddToList(&apply);
     1114
     1115    MTaskEnv taskenv("ExtractSignal");
     1116    taskenv.SetDefault(fExtractor ? fExtractor : &extract2);
     1117
     1118    tlist.AddToList(&taskenv);
     1119    tlist.AddToList(&pinext);
     1120    tlist.AddToList(&blindext);
     1121
     1122    MTaskEnv taskenv2("ExtractTime");
     1123    taskenv2.SetDefault(fTimeExtractor ? fTimeExtractor : &timespline);
     1124
     1125    if (fRelTimes)
     1126        tlist.AddToList(&taskenv2);
     1127
     1128    if (fColor == MCalibrationCam::kCT1)
     1129        tlist.AddToList(&cont);
     1130
     1131    tlist.AddToList(&fillcam);
     1132    tlist.AddToList(&calcalc);
     1133
     1134    MCalibrationRelTimeCalc timecalc;
     1135    if (fRelTimes)
     1136    {
     1137        tlist.AddToList(&filltme);
     1138        tlist.AddToList(&timecalc);
     1139    }
     1140
     1141
     1142    // Create and setup the eventloop
     1143    MEvtLoop evtloop(fName);
     1144    evtloop.SetParList(&plist);
     1145    evtloop.SetDisplay(fDisplay);
     1146    evtloop.SetLogStream(fLog);
     1147    if (fEnv)
     1148        evtloop.ReadEnv(*fEnv);
     1149
     1150    // Execute first analysis
     1151    if (!evtloop.Eventloop())
     1152    {
     1153        *fLog << err << GetDescriptor() << ": Failed." << endl;
     1154        return kFALSE;
     1155    }
     1156
     1157    tlist.PrintStatistics();
     1158
     1159    //
     1160    // The next lines are necessary in order to avoid that
     1161    // the last entry drawn by MFillH gets deleted again from
     1162    // the display. No idea where this comes from...
     1163    //
     1164    /*
     1165    if (fDisplay)
     1166    {
     1167        if (IsUsePINDiode())
    11251168        {
    1126           *fLog << warn << GetDescriptor()
    1127                 << ": No extractor has been chosen, take default MTimeExtractSpline " << endl;
    1128           tlist.AddToList(&timespline);
     1169            MHCalibrationChargePINDiode *pin =
     1170                (MHCalibrationChargePINDiode*)plist.FindObject("MHCalibrationChargePINDiode");
     1171            pin->DrawClone(Form("nonew %s",drawoption.Data()));
    11291172        }
    1130     }
    1131 
    1132   if (fColor == MCalibrationCam::kCT1)
    1133     tlist.AddToList(&cont);
    1134 
    1135   tlist.AddToList(&fillcam);
    1136 
    1137   if (fRelTimes)
    1138     {
    1139       tlist.AddToList(&filltme);
    1140       tlist.AddToList(&timecalc);
    1141     }
    1142  
    1143   if (IsUseBlindPixel())
    1144     tlist.AddToList(&fillbnd);
    1145   if (IsUsePINDiode())
    1146     tlist.AddToList(&fillpin);
    1147 
    1148   tlist.AddToList(&calcalc);
    1149 
    1150   // Create and setup the eventloop
    1151   MEvtLoop evtloop(fName);
    1152   evtloop.SetParList(&plist);
    1153   evtloop.SetDisplay(fDisplay);
    1154   evtloop.SetLogStream(fLog);
    1155  
    1156   // Execute first analysis
    1157   if (!evtloop.Eventloop())
    1158     {
    1159       *fLog << err << GetDescriptor() << ": Failed." << endl;
    1160       return kFALSE;
    1161     }
    1162  
    1163   tlist.PrintStatistics();
    1164 
    1165   //
    1166   // The next lines are necessary in order to avoid that
    1167   // the last entry drawn by MFillH gets deleted again from
    1168   // the display. No idea where this comes from...
    1169   //
    1170   if (fDisplay)
    1171     {
    1172       if (IsUsePINDiode())
     1173        else if (IsUseBlindPixel())
    11731174        {
    1174           MHCalibrationChargePINDiode *pin =
    1175             (MHCalibrationChargePINDiode*)plist.FindObject("MHCalibrationChargePINDiode");
    1176           pin->DrawClone(Form("nonew %s",drawoption.Data()));
     1175            MHCalibrationChargeBlindCam *cam =
     1176                (MHCalibrationChargeBlindCam*)plist.FindObject("MHCalibrationChargeBlindCam");
     1177            cam->DrawClone(Form("nonew %s",drawoption.Data()));
    11771178        }
    1178       else if (IsUseBlindPixel())
     1179        else if (fRelTimes)
    11791180        {
    1180           MHCalibrationChargeBlindCam *cam =
    1181             (MHCalibrationChargeBlindCam*)plist.FindObject("MHCalibrationChargeBlindCam");
    1182           cam->DrawClone(Form("nonew %s",drawoption.Data()));
     1181            MHCalibrationRelTimeCam *cam =
     1182                (MHCalibrationRelTimeCam*)plist.FindObject("MHCalibrationRelTimeCam");
     1183            cam->DrawClone(Form("nonew %s",drawoption.Data()));
    11831184        }
    1184       else if (fRelTimes)
     1185        else
    11851186        {
    1186           MHCalibrationRelTimeCam *cam =
    1187             (MHCalibrationRelTimeCam*)plist.FindObject("MHCalibrationRelTimeCam");
    1188           cam->DrawClone(Form("nonew %s",drawoption.Data()));
     1187            MHCalibrationChargeCam *cam =
     1188                (MHCalibrationChargeCam*)plist.FindObject("MHCalibrationChargeCam");
     1189            cam->DrawClone(Form("nonew %s",drawoption.Data()));
    11891190        }
    1190       else
    1191         {
    1192           MHCalibrationChargeCam *cam =
    1193             (MHCalibrationChargeCam*)plist.FindObject("MHCalibrationChargeCam");
    1194           cam->DrawClone(Form("nonew %s",drawoption.Data()));
    1195         }
    1196     }
    1197  
    1198   DisplayResult(plist);
    1199 
    1200   if (!WriteResult())
    1201     return kFALSE;
    1202  
    1203   *fLog << inf << GetDescriptor() << ": Done." << endl;
    1204  
    1205   return kTRUE;
     1191    }
     1192    */
     1193
     1194    DisplayResult(plist);
     1195
     1196    if (!WriteResult())
     1197        return kFALSE;
     1198
     1199    *fLog << inf << GetDescriptor() << ": Done." << endl;
     1200
     1201    return kTRUE;
    12061202}
    12071203
  • trunk/MagicSoft/Mars/mjobs/MJCalibration.h

    r4604 r4609  
    2121#endif
    2222
     23class TEnv;
    2324class MRunIter;
     25class MSequence;
    2426class MParList;
    2527class MPedestalCam;
    2628class MExtractor;
    2729class MExtractTime;
     30
     31class MExtractBlindPixel;
     32class MHCalibrationChargeBlindCam;
     33
    2834class MJCalibration : public MParContainer
    2935{
    3036private:
     37    static const Int_t gkIFAEBoxInaugurationRun;         // Run number of first IFAE box calibration
     38    static const Int_t gkSecondBlindPixelInstallation;   // Run number upon which second blind pixel was installed
     39    static const Int_t gkThirdBlindPixelInstallation;    // Run number upon which third blind pixel was installed
    3140
    32   static const Int_t gkIFAEBoxInaugurationRun;          // Run number of first IFAE box calibration
    33   static const Int_t gkSecondBlindPixelInstallation;    // Run number upon which second blind pixel was installed
    34   static const Int_t gkThirdBlindPixelInstallation;     // Run number upon which third blind pixel was installed 
     41    TString fOutputPath;                                 // Path to the output files
    3542
    36   TString fOutputPath;                                 // Path to the output files
     43    TEnv           *fEnv;
     44    MRunIter       *fRuns;                               // Calibration files
     45    MSequence      *fSequence;                           // Sequence
     46
     47    MExtractor     *fExtractor;                          // Signal extractor
     48    MExtractTime   *fTimeExtractor;                      // Arrival Time extractor
     49
     50    MBadPixelsCam              fBadPixels;               // Bad Pixels cam, can be set from previous runs
     51    MCalibrationChargeCam      fCalibrationCam;          // Calibration conversion factors FADC2Phe
     52    MCalibrationChargeBlindCam fCalibrationBlindCam;     // Calibration from Blind Pixel(s)
     53    MCalibrationChargePINDiode fCalibrationPINDiode;     // Calibration from PIN Diode
     54    MCalibrationQECam          fQECam;                   // Quantum efficiency, can be set from previous runs
     55    MCalibrationRelTimeCam     fRelTimeCam;              // Calibration constants rel. times
     56
     57    MCalibrationCam::PulserColor_t fColor;               // Colour of the pulsed LEDs
     58
     59    enum  Display_t                                      // Possible Display types
     60    {
     61        kFullDisplay,
     62        kDataCheckDisplay,
     63        kNormalDisplay
     64    };
     65
     66    Display_t fDisplayType;                              // Chosen Display type
     67
     68    enum  Device_t                                       // Possible devices for calibration
     69    {
     70        kUseBlindPixel,
     71        kUsePINDiode
     72    };
     73
     74    Byte_t fDevices;                                     // Bit-field for used devices for calibration
     75
     76    Bool_t fRelTimes;                                    // Flag if relative times have to be calibrated
     77    Bool_t fDataCheck;                                   // Flag if the data check is run on raw data
     78    Bool_t fDebug;
     79
     80    void   DisplayResult(MParList &plist);
     81    Bool_t WriteResult();
     82    void   CheckEnv();
     83
     84    // WORKAROUNDS!!!
     85    Bool_t FindColor();
     86    void   InitBlindPixel(MExtractBlindPixel &blindext,
     87                          MHCalibrationChargeBlindCam &blindcam);
     88
     89public:
     90    MJCalibration(const char *name=NULL, const char *title=NULL);
     91    ~MJCalibration();
    3792 
    38   MRunIter       *fRuns;                               // Calibration files
    39   MExtractor     *fExtractor;                          // Signal extractor
    40   MExtractTime   *fTimeExtractor;                      // Arrival Time extractor
     93    const char* GetOutputFile() const;
     94    void SetEnv(const char *env);
     95
     96    MCalibrationChargeCam  &GetCalibrationCam() { return fCalibrationCam; }
     97    MCalibrationRelTimeCam &GetRelTimeCam()     { return fRelTimeCam;     }
     98    MCalibrationQECam      &GetQECam()          { return fQECam;          }
     99    MBadPixelsCam          &GetBadPixels()      { return fBadPixels;      }
     100
     101    Bool_t IsUseBlindPixel() const;
     102    Bool_t IsUsePINDiode()   const;
     103
     104    void SetBadPixels(const MBadPixelsCam &bad)    { bad.Copy(fBadPixels);   }
     105    void SetExtractor(MExtractor* ext)             { fExtractor = ext; }
     106    void SetTimeExtractor(MExtractTime* ext)       { fTimeExtractor = ext; }
     107    void SetQECam(const MCalibrationQECam &qe)     { qe.Copy(fQECam);        }
     108    void SetColor(const MCalibrationCam::PulserColor_t color) { fColor = color; }
     109
     110    void SetInput(MRunIter *iter) { fRuns=iter; }
     111    void SetSequence(MSequence *seq) { fSequence=seq; }
     112    void SetOutputPath(const char *path=".");
     113
     114    // Displays
     115    void SetFullDisplay()      { fDisplayType = kFullDisplay;      }
     116    void SetDataCheckDisplay() { fDisplayType = kDataCheckDisplay; }
     117    void SetNormalDisplay()    { fDisplayType = kNormalDisplay;    }
     118
     119    // Rel. Time
     120    void SetRelTimeCalibration(const Bool_t b=kTRUE) { fRelTimes = b; }
     121
     122    // Data Check
     123    void SetDataCheck(const Bool_t b=kTRUE) { fDataCheck = b; SetDataCheckDisplay(); }
     124
     125    // Debug
     126    void SetDebug(const Bool_t b=kTRUE) { fDebug = b; }
     127
     128    // Devices
     129    void SetUseBlindPixel(const Bool_t b=kTRUE);
     130    void SetUsePINDiode(const Bool_t b=kTRUE);
    41131 
    42   MBadPixelsCam              fBadPixels;               // Bad Pixels cam, can be set from previous runs
    43   MCalibrationChargeCam      fCalibrationCam;          // Calibration conversion factors FADC2Phe
    44   MCalibrationChargeBlindCam fCalibrationBlindCam;     // Calibration from Blind Pixel(s)
    45   MCalibrationChargePINDiode fCalibrationPINDiode;     // Calibration from PIN Diode
    46   MCalibrationQECam          fQECam;                   // Quantum efficiency, can be set from previous runs
    47   MCalibrationRelTimeCam     fRelTimeCam;              // Calibration constants rel. times
     132    Bool_t ReadCalibrationCam();
     133    Bool_t ProcessFile(MPedestalCam &pedcam);
     134    Bool_t Process(MPedestalCam &pedcam);
    48135
    49   MCalibrationCam::PulserColor_t fColor;               // Colour of the pulsed LEDs
    50 
    51   enum  Display_t   { kFullDisplay, kDataCheckDisplay, kNormalDisplay }; // Possible Display types
    52  
    53   Display_t fDisplayType;                              // Chosen Display type
    54 
    55   enum  Device_t    { kUseBlindPixel, kUsePINDiode  }; // Possible devices for calibration
    56 
    57   Byte_t fDevices;                                     // Bit-field for used devices for calibration
    58  
    59   Bool_t fRelTimes;                                    // Flag if relative times have to be calibrated
    60   Bool_t fDataCheck;                                   // Flag if the data check is run on raw data
    61   Bool_t fDebug;
    62 
    63   void   DisplayResult(MParList &plist);
    64   Bool_t WriteResult();
    65   Bool_t FindColor();
    66  
    67 public:
    68 
    69   MJCalibration(const char *name=NULL, const char *title=NULL);
    70   ~MJCalibration() {}
    71  
    72   const char* GetOutputFile() const;
    73 
    74   MCalibrationChargeCam  &GetCalibrationCam()     { return fCalibrationCam; } 
    75   MCalibrationRelTimeCam &GetRelTimeCam()         { return fRelTimeCam;     }
    76   MCalibrationQECam      &GetQECam()              { return fQECam;          }   
    77   MBadPixelsCam          &GetBadPixels()          { return fBadPixels;      }
    78 
    79   Bool_t IsUseBlindPixel() const;
    80   Bool_t IsUsePINDiode()   const;
    81  
    82   void SetBadPixels(const MBadPixelsCam &bad)     { bad.Copy(fBadPixels);   }
    83   void SetExtractor(MExtractor* ext)              { fExtractor = ext; }
    84   void SetTimeExtractor(MExtractTime* ext)       { fTimeExtractor = ext; }
    85   void SetQECam    (const MCalibrationQECam &qe) { qe.Copy(fQECam);        }   
    86   void SetColor    (const MCalibrationCam::PulserColor_t color) { fColor = color; }
    87 
    88   void SetInput(MRunIter *iter) { fRuns=iter; }
    89   void SetOutputPath(const char *path=".");
    90  
    91   // Displays
    92   void SetFullDisplay()      { fDisplayType = kFullDisplay;      }
    93   void SetDataCheckDisplay() { fDisplayType = kDataCheckDisplay; }
    94   void SetNormalDisplay()    { fDisplayType = kNormalDisplay;    }
    95  
    96   // Rel. Time
    97   void SetRelTimeCalibration(const Bool_t b=kTRUE) { fRelTimes         = b; }
    98 
    99   // Data Check
    100   void SetDataCheck         (const Bool_t b=kTRUE) { fDataCheck        = b; SetDataCheckDisplay(); }
    101 
    102   // Debug
    103   void SetDebug             (const Bool_t b=kTRUE) { fDebug           = b; }
    104  
    105   // Devices
    106   void SetUseBlindPixel( const Bool_t b=kTRUE );
    107   void SetUsePINDiode  ( const Bool_t b=kTRUE ); 
    108  
    109   Bool_t ReadCalibrationCam();
    110   Bool_t ProcessFile( MPedestalCam &pedcam );
    111   Bool_t Process    ( MPedestalCam &pedcam );
    112  
    113   ClassDef(MJCalibration, 0) // Tool to run a calibration per pulser colour and intensity
     136    ClassDef(MJCalibration, 0) // Tool to run a calibration per pulser colour and intensity
    114137};
    115138
  • trunk/MagicSoft/Mars/mjobs/MJPedestal.cc

    r4601 r4609  
    104104{
    105105    if (fSequence)
    106         return Form("%s/calped%05d", (const char*)fOutputPath, fSequence->GetSequence());
     106        return Form("%s/calped%06d.root", (const char*)fOutputPath, fSequence->GetSequence());
    107107
    108108    if (!fRuns)
     
    422422
    423423    return kTRUE;
     424}
     425
     426void MJPedestal::CheckEnv()
     427{
     428    if (!fEnv)
     429        return;
     430
     431    TString e1 = fEnv->GetValue("MJPedestal.OutputPath", "");
     432    if (!e1.IsNull())
     433    {
     434        e1.ReplaceAll("\015", "");
     435        SetOutputPath(e1);
     436    }
    424437}
    425438
     
    449462    *fLog << endl;
    450463
     464    CheckEnv();
     465
    451466    MParList  plist;
    452467    MTaskList tlist;
     
    493508        pedcalc.SetRange(fExtractor->GetHiGainFirst(), fExtractor->GetHiGainLast());
    494509    }
     510  /*
    495511    else
    496512    {
     
    498514        *fLog << ": No extractor has been chosen, take default number of FADC samples " << endl;
    499515    }
     516    */
    500517
    501518    tlist.AddToList(&geomapl);
  • trunk/MagicSoft/Mars/mjobs/MJPedestal.h

    r4601 r4609  
    5858    void   DisplayOutliers(TH1D *hist) const;
    5959
     60    void   CheckEnv();
     61
    6062public:
    6163    MJPedestal(const char *name=NULL, const char *title=NULL);
  • trunk/MagicSoft/Mars/mjobs/MSequence.h

    r4601 r4609  
    5050    // Getter
    5151    UInt_t GetSequence() const { return fSequence; }
     52    UInt_t GetLastRun() const  { return fLastRun; }
    5253
    5354    ClassDef(MSequence, 0)
  • trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.cc

    r4601 r4609  
    185185      fGeom(NULL), fPedContainerName("MPedestalCam")
    186186{
    187   fName  = name  ? name  : "MPedCalcFromLoGain";
    188   fTitle = title ? title : "Task to calculate pedestals from pedestal runs raw data";
    189  
    190   AddToBranchList("fHiGainPixId");
    191   AddToBranchList("fHiGainFadcSamples");
    192  
    193   SetRange(fgHiGainFirst, fgHiGainLast, fgLoGainFirst, fgLoGainLast);
    194 
    195   SetMaxHiGainVar(fgMaxHiGainVar);
    196   SetPedestalUpdate(kTRUE);
    197   Clear();
     187    fName  = name  ? name  : "MPedCalcFromLoGain";
     188    fTitle = title ? title : "Task to calculate pedestals from pedestal runs raw data";
     189
     190    AddToBranchList("fHiGainPixId");
     191    AddToBranchList("fLoGainPixId");
     192    AddToBranchList("fHiGainFadcSamples");
     193    AddToBranchList("fLoGainFadcSamples");
     194
     195    SetRange();
     196    SetMaxHiGainVar();
     197    SetPedestalUpdate(kTRUE);
     198
     199    Clear();
    198200}
    199201
     
    208210void MPedCalcFromLoGain::Clear(const Option_t *o)
    209211{
    210   fRawEvt    = NULL;
    211   fRunHeader = NULL;
    212   fPedestals = NULL;
    213 
    214   // If the size is yet set, set the size
    215   if (fSumx.GetSize()>0)
    216   {
    217       // Reset contents of arrays.
    218       fSumx.Reset();
    219       fSumx2.Reset();
    220       fSumAB0.Reset();
    221       fSumAB1.Reset();
    222       fNumEventsUsed.Reset();
    223       fTotalCounter.Reset();
    224   }
     212    fRawEvt    = NULL;
     213    fRunHeader = NULL;
     214    fPedestals = NULL;
     215
     216    // If the size is yet set, set the size
     217    if (fSumx.GetSize()>0)
     218    {
     219        // Reset contents of arrays.
     220        fSumx.Reset();
     221        fSumx2.Reset();
     222        fSumAB0.Reset();
     223        fSumAB1.Reset();
     224        fNumEventsUsed.Reset();
     225        fTotalCounter.Reset();
     226    }
    225227}
    226228
     
    235237void MPedCalcFromLoGain::SetRange(Byte_t hifirst, Byte_t hilast, Byte_t lofirst, Byte_t lolast)
    236238{
    237   MExtractor::SetRange(hifirst,hilast,lofirst,lolast);
    238 
    239   //
    240   // Redo the checks if the window is still inside the ranges
    241   //
    242   SetWindowSize(fWindowSizeHiGain,fWindowSizeLoGain);
     239    MExtractor::SetRange(hifirst, hilast, lofirst, lolast);
     240
     241    //
     242    // Redo the checks if the window is still inside the ranges
     243    //
     244    SetWindowSize(fWindowSizeHiGain,fWindowSizeLoGain);
    243245}
    244246
     
    247249void MPedCalcFromLoGain::SetMaxHiGainVar(Byte_t maxvar)
    248250{
    249   fMaxHiGainVar = maxvar;
     251    fMaxHiGainVar = maxvar;
    250252}
    251253
     
    265267void MPedCalcFromLoGain::SetWindowSize(Byte_t windowh, Byte_t windowl)
    266268{
    267   fWindowSizeHiGain = windowh & ~1;
    268   fWindowSizeLoGain = windowl & ~1;
    269 
    270   if (fWindowSizeHiGain != windowh)
    271     *fLog << warn << GetDescriptor() << ": Hi Gain window size has to be even, set to: "
    272           << int(fWindowSizeHiGain) << " samples " << endl;
     269    fWindowSizeHiGain = windowh & ~1;
     270    fWindowSizeLoGain = windowl & ~1;
    273271 
    274   if (fWindowSizeLoGain != windowl)
    275     *fLog << warn << GetDescriptor() << ": Lo Gain window size has to be even, set to: "
    276           << int(fWindowSizeLoGain) << " samples " << endl;
     272    if (fWindowSizeHiGain != windowh)
     273    {
     274        *fLog << warn;
     275        *fLog << GetDescriptor() << ": HiGain window has to be even, set to: ";
     276        *fLog << int(fWindowSizeHiGain) << " samples " << endl;
     277    }
    277278   
    278   const Byte_t availhirange = (fHiGainLast-fHiGainFirst+1) & ~1;
    279   const Byte_t availlorange = (fLoGainLast-fLoGainFirst+1) & ~1;
    280 
    281   if (fWindowSizeHiGain > availhirange)
    282     {
    283       *fLog << warn << GetDescriptor()
    284             << Form(": Hi Gain window size: %2i is bigger than available range: [%2i,%2i]",
    285                     (int)fWindowSizeHiGain, (int)fHiGainFirst, (int)fHiGainLast) << endl;
    286       *fLog << warn << GetDescriptor()
    287             << ": Will set window size to: " << (int)availhirange << endl;
    288       fWindowSizeHiGain = availhirange;
    289     }
     279    if (fWindowSizeLoGain != windowl)
     280    {
     281        *fLog << warn;
     282        *fLog << GetDescriptor() << ": Lo Gain window has to be even, set to: ";
     283        *fLog << int(fWindowSizeLoGain) << " samples " << endl;
     284    }
     285     
     286    const Byte_t availhirange = (fHiGainLast-fHiGainFirst+1) & ~1;
     287    const Byte_t availlorange = (fLoGainLast-fLoGainFirst+1) & ~1;
    290288 
    291   if (fWindowSizeLoGain > availlorange)
    292     {
    293       *fLog << warn << GetDescriptor()
    294             << Form(": Lo Gain window size: %2i is bigger than available range: [%2i,%2i]",
    295                     (int)fWindowSizeLoGain, (int)fLoGainFirst, (int)fLoGainLast) << endl;
    296       *fLog << warn << GetDescriptor()
    297             << ": Will set window size to: " << (int)availlorange << endl;
    298       fWindowSizeLoGain = availlorange;
    299     }
    300  
    301   fNumHiGainSamples = (Float_t)fWindowSizeHiGain;
    302   fNumLoGainSamples = (Float_t)fWindowSizeLoGain;
    303  
    304   fSqrtHiGainSamples = TMath::Sqrt(fNumHiGainSamples);
    305   fSqrtLoGainSamples = TMath::Sqrt(fNumLoGainSamples);
     289    if (fWindowSizeHiGain > availhirange)
     290    {
     291        *fLog << warn;
     292        *fLog << GetDescriptor() << ": HiGain window " << (int)fWindowSizeHiGain;
     293        *fLog << " out of range [" << (int)fHiGainFirst;
     294        *fLog << "," << (int)fHiGainLast << "]" << endl;
     295        *fLog << "Will set window size to " << (int)availhirange << endl;
     296        fWindowSizeHiGain = availhirange;
     297    }
     298   
     299    if (fWindowSizeLoGain > availlorange)
     300    {
     301        *fLog << warn;
     302        *fLog << GetDescriptor() << ": LoGain window " << (int)fWindowSizeLoGain;
     303        *fLog << " out of range [" << (int)fLoGainFirst;
     304        *fLog << "," << (int)fLoGainLast << "]" << endl;
     305        *fLog << "Will set window size to " << (int)availlorange << endl;
     306        fWindowSizeLoGain = availlorange;
     307    }
     308    /*
     309     fNumHiGainSamples = (Float_t)fWindowSizeHiGain;
     310     fNumLoGainSamples = (Float_t)fWindowSizeLoGain;
     311
     312     fSqrtHiGainSamples = TMath::Sqrt(fNumHiGainSamples);
     313     fSqrtLoGainSamples = TMath::Sqrt(fNumLoGainSamples);
     314    */
    306315}
    307316
     
    369378Bool_t MPedCalcFromLoGain::ReInit(MParList *pList)
    370379{
    371   Int_t lastdesired   = (Int_t)fLoGainLast;
    372   Int_t lastavailable = (Int_t)fRunHeader->GetNumSamplesLoGain()-1;
    373  
    374   if (lastdesired > lastavailable)
    375     {
    376       const Int_t diff = lastdesired - lastavailable;
    377       *fLog << endl;
    378       *fLog << warn << GetDescriptor()
    379             << Form(": Selected Lo Gain FADC Window [%2i,%2i] ranges out of the available limits: [0,%2i].",
    380                     (int)fLoGainFirst, lastdesired, lastavailable) << endl;
    381       *fLog << GetDescriptor() << ": Will reduce the upper edge to " << (int)(fLoGainLast - diff) << endl;
    382       SetRange(fHiGainFirst, fHiGainLast, fLoGainFirst, fLoGainLast-diff);
    383     }
    384 
    385   lastdesired   = (Int_t)fHiGainLast;
    386   lastavailable = (Int_t)fRunHeader->GetNumSamplesHiGain()-1;
    387  
    388   if (lastdesired > lastavailable)
    389     {
    390       const Int_t diff = lastdesired - lastavailable;
    391       *fLog << endl;
    392       *fLog << warn << GetDescriptor()
    393             << Form(": Selected Hi Gain Range [%2i,%2i] ranges out of the available limits: [0,%2i].",
    394                     (int)fHiGainFirst, lastdesired, lastavailable) << endl;
    395       *fLog << warn << GetDescriptor()
    396             << Form(": Will possibly use %2i samples from the Low-Gain for the High-Gain range", diff)
    397             << endl;
    398       fHiGainLast -= diff;
    399       fHiLoLast    = diff;
    400     }
    401 
    402   lastdesired   = (Int_t)fHiGainFirst+fWindowSizeHiGain-1;
    403   lastavailable = (Int_t)fRunHeader->GetNumSamplesHiGain()-1;
    404  
    405   if (lastdesired > lastavailable)
    406     {
    407       const Int_t diff = lastdesired - lastavailable;
    408       *fLog << endl;
    409       *fLog << warn << GetDescriptor()
    410             << Form(": Selected Hi Gain FADC Window size %2i ranges out of the available limits: [0,%2i].",
    411                     (int)fWindowSizeHiGain, lastavailable) << endl;
    412       *fLog << warn << GetDescriptor()
    413             << Form(": Will use %2i samples from the Low-Gain for the High-Gain extraction", diff)
    414             << endl;
    415 
    416       if ((Int_t)fWindowSizeHiGain > diff)
    417         {
    418           fWindowSizeHiGain -= diff;
    419           fWindowSizeLoGain += diff;
    420         }
    421       else
    422         {
    423           fWindowSizeLoGain += fWindowSizeHiGain;
    424           fLoGainFirst       = diff-fWindowSizeHiGain;
    425           fWindowSizeHiGain  = 0;
    426         }
    427     }
    428 
    429 
    430   // If the size is not yet set, set the size
    431   if (fSumx.GetSize()==0)
    432   {
    433       const Int_t npixels = fPedestals->GetSize();
    434 
    435       fSumx. Set(npixels);
    436       fSumx2.Set(npixels);
    437       fSumAB0.Set(npixels);
    438       fSumAB1.Set(npixels);
    439       fNumEventsUsed.Set(npixels);
    440       fTotalCounter.Set(npixels);
    441 
    442       // Reset contents of arrays.
    443       fSumx.Reset();
    444       fSumx2.Reset();
    445       fSumAB0.Reset();
    446       fSumAB1.Reset();
    447       fNumEventsUsed.Reset();
    448       fTotalCounter.Reset();
    449   }
    450  
    451   if (fWindowSizeHiGain==0 && fWindowSizeLoGain==0)
    452   {
    453       *fLog << err << GetDescriptor()
    454             << ": Number of extracted Slices is 0, abort ... " << endl;
    455       return kFALSE;
    456   }
    457 
    458   *fLog << endl;
    459   *fLog << inf << GetDescriptor() << ": Taking " << (int)fWindowSizeHiGain
    460         << " HiGain FADC samples starting with slice: " << (int)fHiGainFirst << endl;
    461   *fLog << inf << GetDescriptor() << ": Taking " << (int)fWindowSizeLoGain
    462         << " LoGain FADC samples starting with slice: " << (int)fLoGainFirst << endl;
    463  
    464   return kTRUE;
     380    if (fRunHeader->GetNumSamplesHiGain()<1)
     381    {
     382        *fLog << err << "ERROR - Number of available HiGainSamples<1... abort." << endl;
     383        return kFALSE;
     384    }
     385    if (fRunHeader->GetNumSamplesLoGain()<1)
     386    {
     387        *fLog << err << "ERROR - Number of available LoGainSamples<1... abort." << endl;
     388        return kFALSE;
     389    }
     390
     391    fHiGainFirst = TMath::Max(fHiGainFirst, 0);
     392    fHiGainLast  = TMath::Min(fHiGainLast,  fRunHeader->GetNumSamplesHiGain()-1);
     393
     394    fLoGainFirst = TMath::Max(fLoGainFirst, 0);
     395    fLoGainLast  = TMath::Min(fLoGainLast,  fRunHeader->GetNumSamplesLoGain()-1);
     396
     397    fWindowSizeHiGain = TMath::Min(fWindowSizeHiGain, fHiGainLast-fHiGainFirst+1);
     398    fWindowSizeLoGain = TMath::Min(fWindowSizeLoGain, fLoGainLast-fLoGainFirst+1);
     399
     400    SetRange(fHiGainFirst, fHiGainLast, fLoGainFirst, fLoGainLast);
     401
     402    *fLog << inf << endl;
     403    *fLog << "Taking " << Form("%2d", (int)fWindowSizeHiGain) << " HiGain from " << (int)fHiGainFirst << endl;
     404    *fLog << "Taking " << Form("%2d", (int)fWindowSizeLoGain) << " LoGain from " << (int)fLoGainFirst << endl;
     405
     406    if (fWindowSizeHiGain==0)
     407    {
     408        *fLog << err << "ERROR - HiGain windows size == 0... abort." << endl;
     409        return kFALSE;
     410    }
     411    if (fWindowSizeLoGain==0)
     412    {
     413        *fLog << err << "ERROR - HiGain windows size == 0... abort." << endl;
     414        return kFALSE;
     415    }
     416
     417    // If the size is not yet set, set the size
     418    if (fSumx.GetSize()==0)
     419    {
     420        const Int_t npixels = fPedestals->GetSize();
     421
     422        fSumx. Set(npixels);
     423        fSumx2.Set(npixels);
     424        fSumAB0.Set(npixels);
     425        fSumAB1.Set(npixels);
     426        fNumEventsUsed.Set(npixels);
     427        fTotalCounter.Set(npixels);
     428
     429        // Reset contents of arrays.
     430        fSumx.Reset();
     431        fSumx2.Reset();
     432        fSumAB0.Reset();
     433        fSumAB1.Reset();
     434        fNumEventsUsed.Reset();
     435        fTotalCounter.Reset();
     436    }
     437
     438    return kTRUE;
    465439}
    466440
     
    515489        // Find the maximum and minimum signal per slice in the high gain window
    516490        do {
    517             if (*ptr > max) {
     491            if (*ptr > max)
    518492                max = *ptr;
    519             }
    520             if (*ptr < min) {
     493            if (*ptr < min)
    521494                min = *ptr;
    522             }
    523495        } while (++ptr != end);
    524496
     
    527499            continue;
    528500
    529         Byte_t *firstSlice = pixel.GetLoGainSamples() + fLoGainFirst;
    530         Byte_t *lastSlice  = firstSlice + fWindowSizeLoGain;
    531 
    532         Byte_t *slice = firstSlice;
     501        ptr = pixel.GetLoGainSamples() + fLoGainFirst;
     502        end = ptr + fWindowSizeLoGain;
     503
     504        Byte_t *firstSlice = ptr;
     505
    533506        do {
    534             sum += *slice;
    535             sqr += *slice * *slice;
    536         } while (++slice != lastSlice);
     507            sum += *ptr;
     508            sqr += *ptr * *ptr;
     509        } while (++ptr != end);
    537510
    538511        const Float_t msum   = (Float_t)sum;
     
    569542
    570543    if (fPedestalUpdate)
     544    {
     545        fPedestals->ReCalc(*fGeom);
    571546        fPedestals->SetReadyToSave();
     547    }
    572548
    573549    return kTRUE;
     
    581557{
    582558    // Compute pedestals and rms from the whole run
    583     if (fPedestalUpdate)
     559    if (fPedestalUpdate || GetNumExecutions()<1)
    584560        return kTRUE;
    585561
    586     *fLog << flush << inf << "Calculating Pedestals..." << flush;
     562    *fLog << flush << inf << "Calculating pedestals..." << flush;
     563
     564    Double_t sum = 0;
    587565
    588566    const Int_t npix = fGeom->GetNumPixels();
     
    592570        if (n>1)
    593571            Calc(n, idx);
    594     }
    595 
     572        sum += n;
     573    }
     574
     575    *fLog << flush << inf << "Calculating means..." << flush;
     576
     577    fPedestals->SetTotalEntries((UInt_t)(sum/npix*(fWindowSizeLoGain+fWindowSizeHiGain)));
     578    fPedestals->ReCalc(*fGeom);
    596579    fPedestals->SetReadyToSave();
    597580    return kTRUE;
     
    621604        SetWindowSize(hw, lw);
    622605
    623     Int_t num = fNumEventsDump;
    624606    if (IsEnvDefined(env, prefix, "NumEventsDump", print))
    625607    {
    626         num = GetEnvValue(env, prefix, "NumEventsDump", num);
     608        SetDumpEvents(GetEnvValue(env, prefix, "NumEventsDump", fNumEventsDump));
    627609        rc = kTRUE;
    628610    }
    629     SetDumpEvents(num);
    630 
    631     Byte_t max = fMaxHiGainVar;
     611
    632612    if (IsEnvDefined(env, prefix, "MaxHiGainVar", print))
    633613    {
    634         max = GetEnvValue(env, prefix, "MaxHiGainVar", max);
     614        SetMaxHiGainVar(GetEnvValue(env, prefix, "MaxHiGainVar", fMaxHiGainVar));
    635615        rc = kTRUE;
    636616    }
    637     SetMaxHiGainVar(max);
    638 
    639     Bool_t upd = fPedestalUpdate;
     617
    640618    if (IsEnvDefined(env, prefix, "PedestalUpdate", print))
    641619    {
    642         upd = GetEnvValue(env, prefix, "PedestalUpdate", upd);
     620        SetPedestalUpdate(GetEnvValue(env, prefix, "PedestalUpdate", fPedestalUpdate));
    643621        rc = kTRUE;
    644622    }
    645     SetPedestalUpdate(upd);
    646623
    647624    return rc;
  • trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.h

    r4601 r4609  
    6262
    6363    // Setter
    64     void SetRange(Byte_t hifirst=0, Byte_t hilast=0, Byte_t lofirst=0, Byte_t lolast=0);
    65     void SetWindowSize(Byte_t windowh=0, Byte_t windowl=0);
    66     void SetMaxHiGainVar(Byte_t maxvar=0);
     64    void SetRange(Byte_t hifirst=fgHiGainFirst, Byte_t hilast=fgHiGainLast, Byte_t lofirst=fgLoGainFirst, Byte_t lolast=fgLoGainLast);
     65    void SetWindowSize(Byte_t windowh=fgHiGainWindowSize, Byte_t windowl=fgLoGainWindowSize);
     66    void SetMaxHiGainVar(Byte_t maxvar=fgMaxHiGainVar);
    6767    void SetDumpEvents(UInt_t dumpevents = 0) {fNumEventsDump = dumpevents;}
    6868    void SetPedestalUpdate(Bool_t pedupdate)  {fPedestalUpdate = pedupdate;}
  • trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.cc

    r4601 r4609  
    394394
    395395  const Int_t npixels  = fPedestals->GetSize();
    396   const Int_t areas    = fPedestals->GetAverageAreas();
    397   const Int_t sectors  = fPedestals->GetAverageSectors();
     396//  const Int_t areas    = fPedestals->GetAverageAreas();
     397//  const Int_t sectors  = fPedestals->GetAverageSectors();
    398398 
    399399  if (fSumx.GetSize()==0)
     
    402402      fSumx2.Set(npixels);
    403403     
    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);
     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);
    411411     
    412412      fSumx.Reset();
     
    439439Int_t MPedCalcPedRun::Process()
    440440{
     441    MRawEvtPixelIter pixel(fRawEvt);
     442
     443    while (pixel.Next())
     444    {
     445        const UInt_t idx    = pixel.GetPixelId();
     446
     447        Byte_t *ptr = pixel.GetHiGainSamples() + fHiGainFirst;
     448        Byte_t *end = ptr + fWindowSizeHiGain;
     449
     450        UInt_t sum = 0;
     451        UInt_t sqr = 0;
     452
     453        if (fWindowSizeHiGain != 0)
     454        {
     455            do
     456            {
     457                sum += *ptr;
     458                sqr += *ptr * *ptr;
     459            }
     460            while (++ptr != end);
     461        }
     462
     463        if (fWindowSizeLoGain != 0)
     464        {
     465            ptr = pixel.GetLoGainSamples() + fLoGainFirst;
     466            end = ptr + fWindowSizeLoGain;
     467
     468            do
     469            {
     470                sum += *ptr;
     471                sqr += *ptr * *ptr;
     472            }
     473            while (++ptr != end);
     474        }
     475
     476        const Float_t msum = (Float_t)sum;
     477
     478        //
     479        // These three lines have been uncommented by Markus Gaug
     480        // If anybody needs them, please contact me!!
     481        //
     482        //      const Float_t higainped = msum/fNumHiGainSlices;
     483        //      const Float_t higainrms = TMath::Sqrt((msqr-msum*msum/fNumHiGainSlices)/(fNumHiGainSlices-1.));
     484        //      (*fPedestals)[idx].Set(higainped, higainrms);
     485
     486        fSumx[idx]          += msum;
     487        //
     488        // The old version:
     489        //
     490        //       const Float_t msqr = (Float_t)sqr;
     491        //      fSumx2[idx] += msqr;
     492        //
     493        // The new version:
     494        //
     495        const Float_t sqrsum  = msum*msum;
     496        fSumx2[idx]          += sqrsum;
     497    }
     498
     499    fNumSamplesTot += fWindowSizeHiGain + fWindowSizeLoGain;
     500
     501    return kTRUE;
     502/*
    441503  MRawEvtPixelIter pixel(fRawEvt);
    442504 
     
    445507      const UInt_t idx    = pixel.GetPixelId();
    446508      const UInt_t aidx   = (*fGeom)[idx].GetAidx();
    447       const UInt_t sector = (*fGeom)[idx].GetSector();     
     509      const UInt_t sector = (*fGeom)[idx].GetSector();
    448510
    449511      Byte_t *ptr = pixel.GetHiGainSamples() + fHiGainFirst;
     
    490552      fSumx[idx]          += msum;
    491553      fAreaSumx[aidx]     += msum;
    492       fSectorSumx[sector] += msum;     
     554      fSectorSumx[sector] += msum;
    493555      //
    494556      // The old version:
     
    502564      fSumx2[idx]          += sqrsum;
    503565      fAreaSumx2[aidx]     += sqrsum;
    504       fSectorSumx2[sector] += sqrsum;     
     566      fSectorSumx2[sector] += sqrsum;
    505567    }
    506568 
    507569  fPedestals->SetReadyToSave();
    508570  fNumSamplesTot += fWindowSizeHiGain + fWindowSizeLoGain;
    509 
     571*/
    510572  return kTRUE;
    511573}
     
    517579Int_t MPedCalcPedRun::PostProcess()
    518580{
    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     }
     581    const ULong_t nevts = GetNumExecutions();
     582    if (nevts<1)
     583    {
     584        *fLog << err << "ERROR - Not enough events recorded...abort." << endl;
     585        return kFALSE;
     586    }
     587
     588    *fLog << flush << inf << "Calculating pedestals..." << flush;
     589
     590    // Compute pedestals and rms from the whole run
     591    const ULong_t n     = fNumSamplesTot;
     592    const ULong_t npix  = fGeom->GetNumPixels();
     593
     594    for (UInt_t pixidx=0; pixidx<npix; pixidx++)
     595    {
     596        const Float_t sum  = fSumx.At(pixidx);
     597        const Float_t sum2 = fSumx2.At(pixidx);
     598        const Float_t higainped = sum/n;
     599        //
     600        // The old version:
     601        //
     602        //      const Float_t higainrms = TMath::Sqrt((sum2-sum*sum/n)/(n-1.));
     603        //
     604        // The new version:
     605        //
     606        // 1. Calculate the Variance of the sums:
     607        Float_t higainVar = (sum2-sum*sum/nevts)/(nevts-1.);
     608        // 2. Scale the variance to the number of slices:
     609        higainVar /= (Float_t)(fWindowSizeHiGain+fWindowSizeLoGain);
     610        // 3. Calculate the RMS from the Variance:
     611        (*fPedestals)[pixidx].Set(higainped, higainVar < 0 ? 0. : TMath::Sqrt(higainVar), 0, nevts);
     612    }
     613
     614    *fLog << flush << inf << "Calculating means..." << flush;
     615
     616    fPedestals->SetTotalEntries(fNumSamplesTot);
     617    fPedestals->ReCalc(*fGeom);
     618    fPedestals->SetReadyToSave();
     619
     620    return kTRUE;
    553621
    554622  //
    555623  // Loop over the (two) area indices to get the averaged pedestal per aidx
    556624  //
     625  /*
     626
     627   // Compute pedestals and rms from the whole run
     628    const ULong_t n     = fNumSamplesTot;
     629    const ULong_t nevts = GetNumExecutions();
     630
     631    MRawEvtPixelIter pixel(fRawEvt);
     632
     633    while (pixel.Next())
     634    {
     635        const Int_t  pixid  = pixel.GetPixelId();
     636        const UInt_t aidx   = (*fGeom)[pixid].GetAidx();
     637        const UInt_t sector = (*fGeom)[pixid].GetSector();
     638
     639        fAreaValid  [aidx]++;
     640        fSectorValid[sector]++;
     641
     642        const Float_t sum  = fSumx.At(pixid);
     643        const Float_t sum2 = fSumx2.At(pixid);
     644        const Float_t higainped = sum/n;
     645        //
     646        // The old version:
     647        //
     648        //      const Float_t higainrms = TMath::Sqrt((sum2-sum*sum/n)/(n-1.));
     649        //
     650        // The new version:
     651        //
     652        // 1. Calculate the Variance of the sums:
     653        Float_t higainVar = (sum2-sum*sum/nevts)/(nevts-1.);
     654        // 2. Scale the variance to the number of slices:
     655        higainVar /= (Float_t)(fWindowSizeHiGain+fWindowSizeLoGain);
     656        // 3. Calculate the RMS from the Variance:
     657        (*fPedestals)[pixid].Set(higainped, higainVar < 0 ? 0. : TMath::Sqrt(higainVar), 0, nevts);
     658    }
    557659  for (Int_t aidx=0; aidx<fAreaValid.GetSize(); aidx++)
    558660    {
     
    582684      fPedestals->GetAverageArea(aidx).Set(higainped, higainrms);
    583685    }
    584  
     686    */
    585687  //
    586688  // Loop over the (six) sector indices to get the averaged pedestal per sector
    587689  //
     690  /*
    588691  for (Int_t sector=0; sector<fSectorValid.GetSize(); sector++)
    589692    {
     
    613716      fPedestals->GetAverageSector(sector).Set(higainped, higainrms);
    614717    }
    615  
    616718  fPedestals->SetTotalEntries(fNumSamplesTot);
    617719  fPedestals->SetReadyToSave();
    618720
    619721  return kTRUE;
     722    */
    620723}
    621724
  • trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.h

    r4601 r4609  
    3333    TArrayD fSumx;             // sum of values
    3434    TArrayD fSumx2;            // sum of squared values
     35    /*
    3536    TArrayD fAreaSumx;         // averaged sum of values per area idx
    3637    TArrayD fAreaSumx2;        // averaged sum of squared values per area idx
     
    3940    TArrayD fSectorSumx2;      // averaged sum of squared values per sector
    4041    TArrayI fSectorValid;      // number of valid pixel with sector idx
     42    */
    4143
    4244    Int_t  PreProcess (MParList *pList);
  • trunk/MagicSoft/Mars/mpedestal/MPedPhotCam.cc

    r4384 r4609  
    229229void MPedPhotCam::ReCalc(const MGeomCam &geom, MBadPixelsCam *bad)
    230230{
    231    
    232231    const Int_t np = GetSize();
    233232    const Int_t ns = GetNumSectors();
    234233    const Int_t na = GetNumAreas();
    235 
    236    
    237234
    238235    TArrayI acnt(na);
     
    243240    TArrayD ssumr(ns);
    244241
    245 
    246242    for (int i=0; i<np; i++)
    247243    {
    248        
    249 
    250       if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
    251         continue; //was: .IsBad()       
     244        if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
     245            continue; //was: .IsBad()
    252246
    253247        // Create sums for areas and sectors
     
    268262        ssumr[sect] += ne*rms;
    269263        scnt[sect]  += ne;
    270    
    271 
    272     }
    273 
    274     for (int i=0; i<ns; i++){
    275       if (scnt[i]>0)  GetSector(i).Set(ssumx[i]/scnt[i], ssumr[i]/scnt[i], scnt[i]);
    276       else GetSector(i).Set(-1., -1., 0);
    277     }
    278 
    279     for (int i=0; i<na; i++){
    280       if (acnt[i]>0) GetArea(i).Set(asumx[i]/acnt[i], asumr[i]/acnt[i], acnt[i]);
    281       else  GetArea(i).Set(-1., -1., 0);
    282     }
     264    }
     265
     266    for (int i=0; i<ns; i++)
     267        if (scnt[i]>0)
     268            GetSector(i).Set(ssumx[i]/scnt[i], ssumr[i]/scnt[i], scnt[i]);
     269        else
     270            GetSector(i).Clear();
     271
     272    for (int i=0; i<na; i++)
     273        if (acnt[i]>0)
     274            GetArea(i).Set(asumx[i]/acnt[i], asumr[i]/acnt[i], acnt[i]);
     275        else
     276            GetArea(i).Clear();
    283277}
    284278
  • trunk/MagicSoft/Mars/mpedestal/MPedPhotCam.h

    r4384 r4609  
    5151    void Print(Option_t *o="") const;
    5252
    53     void ReCalc(const MGeomCam &geom, MBadPixelsCam *bad);
     53    void ReCalc(const MGeomCam &geom, MBadPixelsCam *bad=NULL);
    5454
    5555    Bool_t GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type=0) const;
  • trunk/MagicSoft/Mars/mpedestal/MPedestalCam.cc

    r4404 r4609  
    3535#include "MPedestalPix.h"
    3636
     37#include <TArrayI.h>
     38#include <TArrayD.h>
     39
    3740#include <TClonesArray.h>
    3841
     
    4346
    4447#include "MGeomCam.h"
     48#include "MGeomPix.h"
     49
     50#include "MBadPixelsCam.h"
     51#include "MBadPixelsPix.h"
    4552
    4653ClassImp(MPedestalCam);
     
    311318}
    312319
    313 Bool_t MPedestalCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const {
    314 
    315   if (GetSize() <= idx)
    316     return kFALSE;
    317 
    318   if (!(*this)[idx].IsValid())
    319     return kFALSE;
    320 
    321   switch (type) {
     320// --------------------------------------------------------------------------
     321//
     322// Calculates the avarage pedestal and pedestal rms for all sectors
     323// and pixel sizes. The geometry container is used to get the necessary
     324// geometry information (sector number, size index) If the bad pixel
     325// container is given all pixels which have the flag 'bad' are ignored
     326// in the calculation of the sector and size average.
     327//
     328void MPedestalCam::ReCalc(const MGeomCam &geom, MBadPixelsCam *bad)
     329{
     330    const Int_t np = GetSize();
     331    const Int_t ns = geom.GetNumSectors();
     332    const Int_t na = geom.GetNumAreas();
     333
     334    TArrayI acnt(na);
     335    TArrayI scnt(ns);
     336    TArrayD asumx(na);
     337    TArrayD ssumx(ns);
     338    TArrayD asumr(na);
     339    TArrayD ssumr(ns);
     340
     341    for (int i=0; i<np; i++)
     342    {
     343        if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
     344            continue; //was: .IsBad()
     345
     346        // Create sums for areas and sectors
     347        const UInt_t aidx = geom[i].GetAidx();
     348        const UInt_t sect = geom[i].GetSector();
     349
     350        const MPedestalPix &pix = (*this)[i];
     351
     352        const UInt_t  ne   = pix.GetNumEvents();
     353        const Float_t mean = pix.GetPedestal();
     354        const Float_t rms  = pix.GetPedestalRms();
     355
     356        asumx[aidx] += ne*mean;
     357        asumr[aidx] += ne*rms;
     358        acnt[aidx]  += ne;
     359
     360        ssumx[sect] += ne*mean;
     361        ssumr[sect] += ne*rms;
     362        scnt[sect]  += ne;
     363    }
     364
     365    for (int i=0; i<ns; i++)
     366        if (scnt[i]>0)
     367            GetAverageSector(i).Set(ssumx[i]/scnt[i], ssumr[i]/scnt[i], 0, scnt[i]);
     368        else
     369            GetAverageSector(i).Clear();
     370
     371    for (int i=0; i<na; i++)
     372        if (acnt[i]>0)
     373            GetAverageArea(i).Set(asumx[i]/acnt[i], asumr[i]/acnt[i], 0, acnt[i]);
     374        else
     375            GetAverageArea(i).Clear();
     376}
     377
     378Bool_t MPedestalCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
     379{
     380    if (GetSize() <= idx)
     381        return kFALSE;
     382
     383    if (!(*this)[idx].IsValid())
     384        return kFALSE;
     385
     386    switch (type)
     387    {
    322388    case 0:
    323       val = (*this)[idx].GetPedestal();
    324       break;
     389        val = (*this)[idx].GetPedestal();
     390        break;
    325391    case 1:
    326       val = fTotalEntries > 0 ?
    327           (*this)[idx].GetPedestalRms()/TMath::Sqrt((Float_t)fTotalEntries)
    328         : (*this)[idx].GetPedestalError();
    329       break;
     392        val = fTotalEntries > 0 ?
     393            (*this)[idx].GetPedestalRms()/TMath::Sqrt((Float_t)fTotalEntries)
     394            : (*this)[idx].GetPedestalError();
     395        break;
    330396    case 2:
    331       val = (*this)[idx].GetPedestalRms();
    332       break;
     397        val = (*this)[idx].GetPedestalRms();
     398        break;
    333399    case 3:
    334       val = fTotalEntries > 0 ?
    335           (*this)[idx].GetPedestalRms()/TMath::Sqrt((Float_t)fTotalEntries)/2.
    336         : (*this)[idx].GetPedestalRmsError();
    337       break;
     400        val = fTotalEntries > 0 ?
     401            (*this)[idx].GetPedestalRms()/TMath::Sqrt((Float_t)fTotalEntries)/2.
     402            : (*this)[idx].GetPedestalRmsError();
     403        break;
    338404    default:
    339       return kFALSE;
     405        return kFALSE;
    340406    }
    341   return kTRUE;
     407    return kTRUE;
    342408}
    343409
    344410void MPedestalCam::DrawPixelContent(Int_t idx) const
    345411{
    346   *fLog << warn << "MPedestalCam::DrawPixelContent - not available." << endl;
    347 }
     412    *fLog << warn << "MPedestalCam::DrawPixelContent - not available." << endl;
     413}
  • trunk/MagicSoft/Mars/mpedestal/MPedestalCam.h

    r3803 r4609  
    1313class MGeomCam;
    1414class MPedestalPix;
     15class MBadPixelsCam;
    1516
    1617class MPedestalCam : public MParContainer, public MCamEvent
     
    5152  void  InitAverageSectors             ( const UInt_t i      );
    5253
     54  void ReCalc(const MGeomCam &geom, MBadPixelsCam *bad=NULL);
     55
    5356  void Print(Option_t *o="") const;
    5457 
  • trunk/MagicSoft/Mars/mpedestal/MPedestalPix.h

    r4556 r4609  
    3737    Float_t GetPedestalRmsError() const { return fNumEvents>0 ? fPedestalRms/TMath::Sqrt((Float_t)fNumEvents/2) : 0; }
    3838
     39    UInt_t  GetNumEvents() const { return fNumEvents; }
     40
    3941    Bool_t IsValid() const;
    4042
  • trunk/MagicSoft/Mars/mraw/MRawEvtHeader.cc

    r4577 r4609  
    7474//
    7575// UInt_t  fNumTrigLvl2
    76 // ------------------ -
     76// --------------------
    7777//
    7878// Number of second level trigger
     
    9090// Type of Trigger.
    9191// This is a Byte (8 bit) to indicated if any of the pixels
    92 // have a non-negligible low gain (1) or not (0)
     92// have a non-negligible low gain (1) or not (0)
     93//
     94// UInt_t fCalibPattern
     95// --------------------
     96//
     97// Bits 1-16: Pulser slot pattern: 16 LEDs slots (see Diploma Thesis
     98// of Michele)
     99//
     100// Bits 17-24: Power of Continous light source: 256 level
     101//
     102// Bits 25-28: Farbe der Continuous Light: red-amber-green/blue-uv
     103//
     104// Bit 29: Calibration trigger On/off
     105// Bit 30: Pedestal trigger on/off
     106// Bit 31: PIN Diode calibration trigger on/off
     107//
     108//
     109// Class Version 2:
     110// ---------------
     111//  - added fCalibPattern
    93112//
    94113/////////////////////////////////////////////////////////////////////////////
     
    379398    return fTrigPattern[0];
    380399}
     400
     401// --------------------------------------------------------------------------
     402//
     403// return pulser slot patter, see class reference
     404//
     405UShort_t MRawEvtHeader::GetPulserSlotPattern() const
     406{
     407    return fTrigPattern[1] & 0xffff;
     408}
     409
     410// --------------------------------------------------------------------------
     411//
     412// return power of continous light, see class reference
     413//
     414Byte_t MRawEvtHeader::GetPowerOfContLight() const
     415{
     416    return (fTrigPattern[1]<<16) & 0xff;
     417}
     418
     419// --------------------------------------------------------------------------
     420//
     421// return color of continous light, see class reference
     422//
     423MRawEvtHeader::CalibCol_t MRawEvtHeader::GetContLedColor() const
     424{
     425    return (CalibCol_t)((fTrigPattern[1]<<24)&0xf);
     426}
  • trunk/MagicSoft/Mars/mraw/MRawEvtHeader.h

    r4577 r4609  
    2424        kTTPedestal    = 1,
    2525        kTTCalibration = 2
     26    };
     27
     28    enum CalibCol_t {
     29        kColRed   = BIT(0),
     30        kColAmber = BIT(1),
     31        kColGreen = BIT(2),
     32        kColBlue  = BIT(3),
     33        kColUV    = BIT(4)
    2634    };
    2735
     
    6068    void FillHeader(UInt_t, Float_t=0);
    6169
    62     UShort_t GetTrigType() const     { return fTrigType; }
    63     UInt_t   GetNumTrigLvl1() const  { return fNumTrigLvl1; }
    64     UInt_t   GetNumTrigLvl2() const  { return fNumTrigLvl2; }
    65     UInt_t   GetDAQEvtNumber() const { return fDAQEvtNumber; }
     70    UShort_t   GetTrigType() const     { return fTrigType; }
     71    UInt_t     GetNumTrigLvl1() const  { return fNumTrigLvl1; }
     72    UInt_t     GetNumTrigLvl2() const  { return fNumTrigLvl2; }
     73    UInt_t     GetDAQEvtNumber() const { return fDAQEvtNumber; }
    6674   
    67     UInt_t GetTriggerID() const;
     75    UInt_t     GetTriggerID() const;
     76
     77    UShort_t   GetPulserSlotPattern() const;
     78    Byte_t     GetPowerOfContLight() const;
     79    CalibCol_t GetContLedColor() const;
    6880
    6981    Int_t ReadEvt(istream& fin, UShort_t ver);
    70     void SkipEvt(istream& fin, UShort_t ver);
     82    void  SkipEvt(istream& fin, UShort_t ver);
    7183
    7284    ClassDef(MRawEvtHeader, 1) // Parameter Conatiner for raw EVENT HEADER
Note: See TracChangeset for help on using the changeset viewer.