Ignore:
Timestamp:
07/13/05 19:06:26 (19 years ago)
Author:
tbretz
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars/mcalib
Files:
16 edited

Legend:

Unmodified
Added
Removed
  • trunk/MagicSoft/Mars/mcalib/CalibLinkDef.h

    r6663 r7188  
    1919#pragma link C++ class MCalibrationIntensityCam+;
    2020#pragma link C++ class MCalibrationIntensityChargeCam+;
     21#pragma link C++ class MCalibrationIntensityConstCam+;
    2122#pragma link C++ class MCalibrationIntensityBlindCam+;
    2223#pragma link C++ class MCalibrationIntensityQECam+;
    2324#pragma link C++ class MCalibrationIntensityRelTimeCam+;
    2425#pragma link C++ class MCalibrationIntensityTestCam+;
     26
    2527#pragma link C++ class MCalibrationCam+;
    2628#pragma link C++ class MCalibrationPix+;
  • trunk/MagicSoft/Mars/mcalib/MCalibConstCam.cc

    r6457 r7188  
    8585// --------------------------------------------------------------------------
    8686//
    87 // Set the size of the camera
    88 //
    89 void MCalibConstCam::InitSize(const UInt_t i)
    90 {
    91     fArray->ExpandCreate(i);
    92 }
    93 
    94 // -------------------------------------------------------------------
    95 //
    96 // Calls TClonesArray::ExpandCreate() for:
    97 // - fAverageAreas
    98 //
    99 void MCalibConstCam::InitAverageAreas(const UInt_t i)
    100 {
    101   fAverageAreas->ExpandCreate(i);
    102 }
    103 
    104 // -------------------------------------------------------------------
    105 //
    106 // Calls TClonesArray::ExpandCreate() for:
    107 // - fAverageSectors
    108 //
    109 void MCalibConstCam::InitAverageSectors(const UInt_t i)
    110 {
    111   fAverageSectors->ExpandCreate(i);
    112 }
     87// Copy 'constructor'
     88//
     89void MCalibConstCam::Copy(TObject &obj) const
     90{
     91
     92  MParContainer::Copy(obj);
     93
     94  MCalibConstCam &cam = (MCalibConstCam&)obj;
     95 
     96  Int_t n = GetSize();
     97 
     98  if (n==0)
     99    return;
     100 
     101  cam.InitSize(n);
     102  for (int i=0; i<n; i++)
     103    (*this)[i].Copy(cam[i]);
     104
     105  n = GetNumAverageArea();
     106  cam.InitAverageAreas(n);
     107  for (int i=0; i<n; i++)
     108    GetAverageArea(i).Copy(cam.GetAverageArea(i));
     109
     110  n = GetNumAverageSector();
     111  cam.InitAverageSectors(n);
     112  for (int i=0; i<n; i++)
     113    GetAverageSector(i).Copy(cam.GetAverageSector(i));
     114}
     115
    113116
    114117// -------------------------------------------------------------------
  • trunk/MagicSoft/Mars/mcalib/MCalibConstCam.h

    r6073 r7188  
    99#endif
    1010
    11 class TClonesArray;
     11#ifndef ROOT_TClonesArray
     12#include <TClonesArray.h>
     13#endif
     14
    1215
    1316class MGeomCam;
     
    2124  TClonesArray *fAverageSectors;  //-> Array of MCalibConstPix, one per camera sector
    2225
     26  Int_t fRunNumber;               // Run number
     27 
    2328public:
    2429
     
    2732 
    2833  void Clear(Option_t *o="");
     34  void Copy(TObject &object) const;
    2935 
    3036  // Getters
    31         MCalibConstPix &GetAverageArea   ( UInt_t i );
    32   const MCalibConstPix &GetAverageArea   ( UInt_t i )            const;
    33   const Int_t           GetNumAverageArea()                      const;
    34         MCalibConstPix &GetAverageSector ( UInt_t i );
    35   const MCalibConstPix &GetAverageSector ( UInt_t i )            const;
    36   const Int_t           GetNumAverageSector()                    const;
    37   Int_t                 GetSize          ()                      const;
     37        MCalibConstPix &GetAverageArea     ( UInt_t i );
     38  const MCalibConstPix &GetAverageArea     ( UInt_t i )  const;
     39  const Int_t           GetNumAverageArea  ()            const;
     40        MCalibConstPix &GetAverageSector   ( UInt_t i );
     41  const MCalibConstPix &GetAverageSector   ( UInt_t i )  const;
     42  const Int_t           GetNumAverageSector()            const;
     43  Int_t                 GetSize            ()            const;
    3844
    39         MCalibConstPix &operator[]     ( Int_t i             );
    40   const MCalibConstPix &operator[]     ( Int_t i             ) const;
     45        MCalibConstPix &operator[]         ( Int_t i  );
     46  const MCalibConstPix &operator[]         ( Int_t i  ) const;
    4147
    4248  void  Init                           ( const MGeomCam &geom);
    43   void  InitSize                       ( const UInt_t i      );
    44   void  InitAverageAreas               ( const UInt_t i      );
    45   void  InitAverageSectors             ( const UInt_t i      );
     49  void  InitSize                       ( const UInt_t i ) { fArray->ExpandCreate(i);          }
     50  void  InitAverageAreas               ( const UInt_t i ) { fAverageAreas->ExpandCreate(i);   }
     51  void  InitAverageSectors             ( const UInt_t i ) { fAverageSectors->ExpandCreate(i); }
    4652
    4753  void Print(Option_t *o="") const;
     54
     55  // Setters
     56  void SetRunNumber( const Int_t n )   {  fRunNumber = n; }
    4857 
    49   Bool_t GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type=0) const;
    50   void  DrawPixelContent(Int_t idx) const;
     58  Bool_t GetPixelContent (Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type=0) const;
     59  void   DrawPixelContent(Int_t idx) const;
    5160
    52   ClassDef(MCalibConstCam, 0)   // Temporary Storage for calibration constants
     61  ClassDef(MCalibConstCam, 1)   // Temporary Storage for calibration constants
    5362};
    5463
  • trunk/MagicSoft/Mars/mcalib/MCalibConstPix.h

    r6073 r7188  
    99{
    1010private:
    11 
    1211    Float_t fCalibConst;        // conversion factor (modified after each interlaced cal. update)
    1312    Float_t fCalibFFactor;      // global F-Factor   (modified after each interlaced cal. update)
    1413
    1514public:
    16  
    1715    MCalibConstPix();
    18    
     16
     17    // TObject
    1918    void Clear(Option_t *o="")
    20       {
    21         fCalibConst   = -1.;
    22         fCalibFFactor = -1.;
    23       }
    24    
     19    {
     20        fCalibConst   = -1.;
     21        fCalibFFactor = -1.;
     22    }
     23
     24    void Copy(TObject &object) const
     25    {
     26        MCalibConstPix &pix =  (MCalibConstPix&)object;
     27        pix.fCalibConst   = fCalibConst;
     28        pix.fCalibFFactor = fCalibFFactor;
     29    }
     30
    2531    // Getters
    2632    Float_t GetCalibConst()   const { return fCalibConst;   }
     
    3137    void SetCalibFFactor( const Float_t f )  { fCalibFFactor = f; }
    3238
    33     ClassDef(MCalibConstPix, 0) // Temporay Storage Calibraion Constant of one pixel
     39    ClassDef(MCalibConstPix, 1) // Temporay Storage Calibraion Constant of one pixel
    3440};
    3541
  • trunk/MagicSoft/Mars/mcalib/MCalibrationBlindCam.cc

    r5128 r7188  
    7171
    7272}
    73 
    74 // --------------------------------------------------------------------------
    75 //
    76 // Copy 'constructor'
    77 //
    78 void MCalibrationBlindCam::Copy(TObject& object) const
    79 {
    80  
    81   MCalibrationBlindCam &calib = (MCalibrationBlindCam&)object;
    82  
    83   MCalibrationCam::Copy(calib);
    84 
    85   /* 
    86   const UInt_t n = GetSize();
    87   if (n != 0)
    88     {
    89       calib.InitSize(n);
    90       for (UInt_t i=0; i<n; i++)
    91         (*this)[i].Copy(calib[i]);
    92     }
    93   */
    94 }
    95 
    96 
    9773
    9874// --------------------------------------------------------------------------
  • trunk/MagicSoft/Mars/mcalib/MCalibrationBlindCam.h

    r5128 r7188  
    1313
    1414public:
    15 
    1615  MCalibrationBlindCam(Int_t nblind=1,const char *name=NULL, const char *title=NULL);
    1716 
    18   void Copy(TObject& object) const;
    19 
    2017  // Inits
    2118  void  Init ( const MGeomCam &geom ) {}
    2219
     20  // Getter
    2321  Bool_t IsFluxInsidePlexiglassAvailable () const;
    2422
     
    2725  Float_t GetFluxInsidePlexiglassRelVar() const;
    2826
     27  // Setter
    2928  void  SetPulserColor( const MCalibrationCam::PulserColor_t col );
    3029 
  • trunk/MagicSoft/Mars/mcalib/MCalibrationChargePix.cc

    r7099 r7188  
    7373//  + fNumSaturated
    7474//
     75// Class Version 4:
     76//  +  Float_t fConversionHiLoSigma;             // Sigma of conversion factor betw. Hi and Lo Gain
     77//  -  Float_t fMeanConvFADC2PheVar;             // Variance conversion factor (F-factor method)
     78//  +  Float_t fMeanConvFADC2PheStatVar;         // Variance conversion factor, only stat. error
     79//  +  Float_t fMeanConvFADC2PheSystVar;         // Variance conversion factor, only syst. error
     80//  -  Float_t fMeanFFactorFADC2PhotVar;         // Variance mean F-Factor photons (F-factor method)
     81//  +  Float_t fMeanFFactorFADC2PhotVar;         // Variance mean F-Factor photons, only stat. error
     82//  -  Float_t fPheFFactorMethod;                // Number Phe's calculated (F-factor method)
     83//  -  Float_t fPheFFactorMethodVar;             // Variance number of Phe's (F-factor method)
     84//  +  Float_t fPheFFactorMethod;                // Number Phe's calculated  with F-factor method)
     85//  +  Float_t fPheFFactorMethodStatVar;         // Variance number of Phe's, only stat. error
     86//  +  Float_t fPheFFactorMethodSystVar;         // Variance number of Phe's, only syst. error
     87//
    7588/////////////////////////////////////////////////////////////////////////////
    7689#include "MCalibrationChargePix.h"
     
    90103const Float_t MCalibrationChargePix::gkFFactorErr          = 0.02;
    91104
    92 const Float_t MCalibrationChargePix::fgConversionHiLo           = 10.;
    93 const Float_t MCalibrationChargePix::fgConversionHiLoErr        = 2.5;
    94 const Float_t MCalibrationChargePix::fgPheFFactorMethodLimit    = 1.;
    95 const Float_t MCalibrationChargePix::fgConvFFactorRelErrLimit   = 0.85;
     105const Float_t MCalibrationChargePix::fgConversionHiLo         = 10.;
     106const Float_t MCalibrationChargePix::fgConversionHiLoErr      = 0.05;
     107const Float_t MCalibrationChargePix::fgConversionHiLoSigma    = 2.5;
     108const Float_t MCalibrationChargePix::fgPheFFactorMethodLimit  = 1.;
     109const Float_t MCalibrationChargePix::fgConvFFactorRelErrLimit = 0.85;
    96110
    97111// --------------------------------------------------------------------------
     
    157171
    158172  fPheFFactorMethod                 =  -1.;
    159   fPheFFactorMethodVar              =  -1.;
     173  fPheFFactorMethodStatVar          =  -1.;
     174  fPheFFactorMethodSystVar          =  -1.;
    160175
    161176  fMeanConvFADC2Phe                 =  -1.;
    162   fMeanConvFADC2PheVar              =  -1.;
     177  fMeanConvFADC2PheStatVar          =  -1.;
     178  fMeanConvFADC2PheSystVar          =  -1.;
    163179  fMeanFFactorFADC2Phot             =  -1.;
    164180  fMeanFFactorFADC2PhotVar          =  -1.; 
     
    168184  MCalibrationPix::Clear();
    169185}
     186
     187// -----------------------------------------------------
     188//
     189void MCalibrationChargePix::Copy(TObject& object) const
     190{
     191  MCalibrationChargePix &pix = (MCalibrationChargePix&)object;
     192 
     193  MCalibrationPix::Copy(pix);
     194
     195  //
     196  // Copy the data members
     197  //
     198  pix.fAbsTimeMean                = fAbsTimeMean              ;
     199  pix.fAbsTimeRms                 = fAbsTimeRms               ;
     200  pix.fCalibFlags                 = fCalibFlags               ;
     201  pix.fConversionHiLo             = fConversionHiLo           ;
     202  pix.fConversionHiLoVar          = fConversionHiLoVar        ;
     203  pix.fConversionHiLoSigma        = fConversionHiLoSigma      ;
     204  pix.fConvFFactorRelVarLimit     = fConvFFactorRelVarLimit   ;
     205  pix.fLoGainPedRmsSquare         = fLoGainPedRmsSquare       ;
     206  pix.fLoGainPedRmsSquareVar      = fLoGainPedRmsSquareVar    ;
     207  pix.fMeanConvFADC2Phe           = fMeanConvFADC2Phe         ;
     208  pix.fMeanConvFADC2PheStatVar    = fMeanConvFADC2PheStatVar  ;
     209  pix.fMeanConvFADC2PheSystVar    = fMeanConvFADC2PheSystVar  ;
     210  pix.fMeanFFactorFADC2Phot       = fMeanFFactorFADC2Phot     ;
     211  pix.fMeanFFactorFADC2PhotVar    = fMeanFFactorFADC2PhotVar  ;
     212  pix.fPed                        = fPed                      ;
     213  pix.fPedVar                     = fPedVar                   ;
     214  pix.fPedRms                     = fPedRms                     ;
     215  pix.fPedRmsVar                  = fPedRmsVar                  ;
     216  pix.fPheFFactorMethod           = fPheFFactorMethod           ;
     217  pix.fPheFFactorMethodStatVar    = fPheFFactorMethodStatVar    ;
     218  pix.fPheFFactorMethodSystVar    = fPheFFactorMethodSystVar    ;
     219  pix.fPheFFactorMethodLimit      = fPheFFactorMethodLimit      ;
     220  pix.fRSigmaSquare               = fRSigmaSquare               ;
     221  pix.fRSigmaSquareVar            = fRSigmaSquareVar            ;
     222  pix.fNumSaturated               = fNumSaturated               ;
     223}                                 
    170224
    171225
     
    185239void MCalibrationChargePix::SetPedestal(const Float_t ped, const Float_t pedrms, const Float_t pederr)
    186240{
    187 
    188241  fPed       = ped;   
    189242  fPedRms    = pedrms;
     
    197250void MCalibrationChargePix::SetPed(const Float_t ped, const Float_t pederr)
    198251{
    199 
    200252  fPed       = ped;   
    201253  fPedVar    = pederr*pederr;
     
    208260void MCalibrationChargePix::SetPedRMS( const Float_t pedrms, const Float_t pedrmserr)
    209261{
    210  
    211262  fPedRms    = pedrms;
    212263  fPedRmsVar = pedrmserr*pedrmserr;
    213  
    214 }
    215 
    216 
    217 // -------------------------------------------------------------------------------
    218 //
    219 // Get the conversion Error Hi-Gain to Low-Gain:
    220 // - If fConversionHiLoVar is smaller than 0 (i.e. has not yet been set), return -1.
     264}
     265
     266
     267// Inline Functions:
     268// -----------------
     269//
     270// GetConversionHiLoErr:
     271//   Get the conversion Error Hi-Gain to Low-Gain:
     272//   If fConversionHiLoVar is smaller than 0 (i.e. has not yet
     273//   been set), return -1.
    221274// 
    222 Float_t MCalibrationChargePix::GetConversionHiLoErr()  const
    223 {
    224   if (fConversionHiLoVar < 0.)
    225     return -1.;
    226 
    227   return TMath::Sqrt(fConversionHiLoVar);
    228 }
     275// GetPedRms(): Get the pedestals RMS: Test bit kHiGainSaturation:
     276//  If yes, return square root of fLoGainPedRmsSquare (if greater than 0,
     277//  otherwise -1.), If no,  return fPedRms
     278//
     279// GetConvertedMean(): Get the Low Gain Mean Charge converted to High Gain
     280//  amplification: Returns fLoGainMean multiplied with fConversionHiLo if
     281//  IsHiGainSaturation(), else return fHiGainMean
     282//
     283// GetConvertedMeanErr(): Get the Error of the converted Low Gain Mean:
     284//  Returns -1 if the variable fLoGainMean or fLoGainMeanVar are smaller than 0.
     285//  Returns the square root of the quadratic sum of the relative variances of
     286//  the fLoGainMean and fConversionHiLo, mulitplied with GetConvertedMean()
     287//  in case of HiGain Saturation,
     288//  else return GetMeanErr()
     289//
     290// GetConvertedSigma(): Get the Low Gain Sigma converted to High Gain
     291//  amplification: Returns fLoGainSigma multiplied with fConversionHiLo
     292//  if IsHiGainSaturation() else return fHiGainSigma
     293//
     294// GetConvertedSigmaErr(): Get the Error of the converted Sigma:
     295//  Returns -1 if the variable fLoGainSigma or fLoGainSigmaVar are smaller than 0.
     296//  if IsHiGainSaturatio()
     297//  returns the square root of the quadratic sum of the relative variances of
     298//  the fLoGainSigma and fConversionHiLo, mulitplied with GetConvertedSigma()
     299//  else returns GetSigmaErr()
     300//
     301// GetConvertedRSigma(): Get the converted reduced Sigma:
     302//  If fRSigmaSquare is smaller than 0 (i.e. has not yet been set), return -1.
     303//  Test bit kHiGainSaturation:
     304//  If yes, return square root of fRSigmaSquare, multiplied with fConversionHiLo,
     305//  If no , return square root of fRSigmaSquare
     306//
     307// GetRSigma(): Get the reduced Sigma:
     308//  If fRSigmaSquare is smaller than 0 (i.e. has not yet been set), return -1.
     309//
     310// GetConvertedRSigmaSquare(): Get the reduced Sigma Square:
     311//  If fRSigmaSquare is smaller than 0 (i.e. has not yet been set), return -1.
     312//  Test bit kHiGainSaturation:
     313//  If yes, return fRSigmaSquare, multiplied with fConversionHiLo^2,
     314//  If no , return fRSigmaSquare
     315//
     316// GetPheFFactorMethodErr(): Get the error on the number of photo-electrons
     317//  (F-Factor Method): If fPheFFactorMethodVar is smaller than 0 (i.e. has
     318//  not yet been set), return -1. Else returns the square root of
     319//  fPheFFactorMethodVar
     320//
     321// GetMeanFFactorFADC2PhotErr(): Get the error on the mean total F-Factor
     322//  of the signal readout (F-Factor Method): If fMeanFFactorFADC2PhotVar
     323//  is smaller than 0 (i.e. has not yet been set), return -1. Else returns
     324//  the square root of fMeanFFactorFADC2PhotVar
     325//
     326// GetPheFFactorMethodRelVar(): Get the relative variance on the number of
     327//  photo-electrons (F-Factor Method): If fPheFFactorMethodVar is smaller
     328//  than 0 (i.e. has not yet been set), return -1. If fPheFFactorMethod
     329//  is 0, return -1. Else returns fPheFFactorMethodVar / fPheFFactorMethod^2
     330//
     331// GetMeanConvFADC2PheErr(): Get the error on the mean conversion factor
     332//  (FFactor  Method): If fMeanConvFADC2PheVar is smaller than 0 (i.e. has
     333//  not yet been set), return -1. Else returns the square root of
     334//  fMeanConvFADC2PheVar
    229335
    230336// --------------------------------------------------------------------------
     
    278384 
    279385
    280 // --------------------------------------------------------------------------
    281 //
    282 // Get the pedestals RMS:
    283 // - Test bit kHiGainSaturation:
    284 //   If yes, return square root of fLoGainPedRmsSquare (if greater than 0, otherwise -1.),
    285 //   If no,  return fPedRms
    286 //
    287 Float_t  MCalibrationChargePix::GetPedRms()  const
    288 {
    289 
    290   if (IsHiGainSaturation())
    291     if (fLoGainPedRmsSquare < 0.)
    292       return -1.;
    293     else
    294       return TMath::Sqrt(fLoGainPedRmsSquare);
    295  
    296   return fPedRms;
    297 }
    298386
    299387// --------------------------------------------------------------------------
     
    321409// --------------------------------------------------------------------------
    322410//
    323 // Get the Low Gain Mean Charge converted to High Gain amplification:
    324 // Returns fLoGainMean multiplied with fConversionHiLo if IsHiGainSaturation(),
    325 //         else return fHiGainMean
    326 //
    327 Float_t MCalibrationChargePix::GetConvertedMean()  const
    328 {
    329 
    330   if (IsHiGainSaturation())
    331     return fLoGainMean * fConversionHiLo;
    332 
    333   return fHiGainMean;
    334 }
    335 
    336 // --------------------------------------------------------------------------
    337 //
    338 // Get the Error of the converted Low Gain Mean:
    339 //
    340 // Returns -1 if the variable fLoGainMean or fLoGainMeanVar are smaller than 0.
    341 //
    342 // Returns the square root of the quadratic sum of the relative variances of
    343 // the fLoGainMean and fConversionHiLo, mulitplied with GetConvertedMean()
    344 // in case of HiGain Saturation,
    345 // else return GetMeanErr()
    346 //
    347 Float_t MCalibrationChargePix::GetConvertedMeanErr()  const
    348 {
    349 
    350   if (IsHiGainSaturation())
    351     {
    352       const Float_t logainrelvar = GetLoGainMeanRelVar();
    353      
    354       if (logainrelvar < 0.)
    355         return -1.;
    356 
    357       return TMath::Sqrt(logainrelvar + GetConversionHiLoRelVar()) * GetConvertedMean();
    358     }
    359   else
    360     return GetMeanErr();
    361 
    362 }
    363 
    364 // --------------------------------------------------------------------------
    365 //
    366 // Get the Low Gain Sigma converted to High Gain amplification:
    367 // Returns fLoGainSigma multiplied with fConversionHiLo if IsHiGainSaturation()
    368 // else return fHiGainSigma
    369 //
    370 Float_t MCalibrationChargePix::GetConvertedSigma()  const
    371 {
    372  
    373   if (IsHiGainSaturation())
    374     return fLoGainSigma * fConversionHiLo;
    375   else
    376     return fHiGainSigma;
    377 }
    378 
    379 // --------------------------------------------------------------------------
    380 //
    381 // Get the Error of the converted Sigma:
    382 //
    383 // Returns -1 if the variable fLoGainSigma or fLoGainSigmaVar are smaller than 0.
    384 //
    385 // if IsHiGainSaturatio()
    386 // returns the square root of the quadratic sum of the relative variances of
    387 // the fLoGainSigma and fConversionHiLo, mulitplied with GetConvertedSigma()
    388 // else returns GetSigmaErr()
    389 //
    390 Float_t MCalibrationChargePix::GetConvertedSigmaErr()  const
    391 {
    392 
    393   if (IsHiGainSaturation())
    394     {
    395       if (fLoGainSigmaVar < 0.)
    396         return -1.;
    397      
    398       if (fLoGainSigma < 0.)
    399         return -1.;
    400      
    401       const Float_t sigmaRelVar =  fLoGainSigmaVar
    402                                 /( fLoGainSigma * fLoGainSigma );
    403 
    404       return TMath::Sqrt(sigmaRelVar+GetConversionHiLoRelVar()) * GetConvertedSigma();
    405     }
    406   else
    407     return GetSigmaErr();
    408 
    409 
    410 }
    411 
    412 // --------------------------------------------------------------------------
    413 //
    414 // Get the converted reduced Sigma:
    415 // - If fRSigmaSquare is smaller than 0 (i.e. has not yet been set), return -1.
    416 // - Test bit kHiGainSaturation:
    417 //   If yes, return square root of fRSigmaSquare, multiplied with fConversionHiLo,
    418 //   If no , return square root of fRSigmaSquare
    419 //
    420 Float_t MCalibrationChargePix::GetConvertedRSigma()  const
    421 {
    422   if (fRSigmaSquare < 0)
    423     return -1;
    424 
    425   const Float_t rsigma = TMath::Sqrt(fRSigmaSquare);
    426  
    427   return IsHiGainSaturation() ? rsigma*fConversionHiLo : rsigma ;
    428 }
    429 
    430 // --------------------------------------------------------------------------
    431 //
    432411// Get the error of the converted reduced Sigma:
    433412// - If fRSigmaSquareVar is smaller than 0 (i.e. has not yet been set), return -1.
     
    460439// --------------------------------------------------------------------------
    461440//
    462 // Get the reduced Sigma:
    463 // - If fRSigmaSquare is smaller than 0 (i.e. has not yet been set), return -1.
    464 //
    465 Float_t MCalibrationChargePix::GetRSigma()  const
    466 {
    467   if (fRSigmaSquare < 0)
    468     return -1;
    469 
    470   return TMath::Sqrt(fRSigmaSquare);
    471  
    472 }
    473 
    474 // --------------------------------------------------------------------------
    475 //
    476441// Get the error of the reduced Sigma:
    477442// - If fRSigmaSquareVar is smaller than 0 (i.e. has not yet been set), return -1.
     
    540505 
    541506  return TMath::Sqrt(rsigmarelvar + meanrelvar) * GetRSigmaPerCharge();
    542 }
    543 
    544 // --------------------------------------------------------------------------
    545 //
    546 // Get the reduced Sigma Square:
    547 // - If fRSigmaSquare is smaller than 0 (i.e. has not yet been set), return -1.
    548 // - Test bit kHiGainSaturation:
    549 //   If yes, return fRSigmaSquare, multiplied with fConversionHiLo^2,
    550 //   If no , return fRSigmaSquare
    551 //
    552 Float_t MCalibrationChargePix::GetConvertedRSigmaSquare()  const
    553 {
    554   if (fRSigmaSquare < 0)
    555     return -1;
    556 
    557   return IsHiGainSaturation() ? fRSigmaSquare*fConversionHiLo*fConversionHiLo : fRSigmaSquare ;
    558507}
    559508
     
    582531}
    583532
    584 // --------------------------------------------------------------------------
    585 //
    586 // Get the error on the number of photo-electrons (F-Factor Method):
    587 // - If fPheFFactorMethodVar is smaller than 0 (i.e. has not yet been set), return -1.
    588 // - Else returns the square root of fPheFFactorMethodVar
    589 //
    590 Float_t MCalibrationChargePix::GetPheFFactorMethodErr()  const
    591 {
    592   if (fPheFFactorMethodVar < 0.)
    593     return -1.;
    594   return TMath::Sqrt(fPheFFactorMethodVar);
    595 }
    596 
    597 // --------------------------------------------------------------------------
    598 //
    599 // Get the error on the mean total F-Factor of the signal readout (F-Factor Method):
    600 // - If fMeanFFactorFADC2PhotVar is smaller than 0 (i.e. has not yet been set), return -1.
    601 // - Else returns the square root of fMeanFFactorFADC2PhotVar
    602 //
    603 Float_t MCalibrationChargePix::GetMeanFFactorFADC2PhotErr()  const
    604 {
    605   if (fMeanFFactorFADC2PhotVar < 0.)
    606     return -1.;
    607   return TMath::Sqrt(fMeanFFactorFADC2PhotVar);
    608 }
    609 
    610 // --------------------------------------------------------------------------
    611 //
    612 // Get the relative variance on the number of photo-electrons (F-Factor Method):
    613 // - If fPheFFactorMethodVar is smaller than 0 (i.e. has not yet been set), return -1.
    614 // - If fPheFFactorMethod    is 0, return -1.
    615 // - Else returns fPheFFactorMethodVar / fPheFFactorMethod^2
    616 //
    617 Float_t MCalibrationChargePix::GetPheFFactorMethodRelVar()  const
    618 {
    619   if (fPheFFactorMethodVar < 0.)
    620     return -1.;
    621   if (fPheFFactorMethod  == 0.)
    622     return -1.;
    623 
    624   return fPheFFactorMethodVar / (fPheFFactorMethod * fPheFFactorMethod);
    625 }
    626 
    627 
    628 // --------------------------------------------------------------------------
    629 //
    630 // Get the error on the mean conversion factor (FFactor  Method):
    631 // - If fMeanConvFADC2PheVar is smaller than 0 (i.e. has not yet been set), return -1.
    632 // - Else returns the square root of fMeanConvFADC2PheVar
    633 //
    634 Float_t MCalibrationChargePix::GetMeanConvFADC2PheErr()  const
    635 {
    636   if (fMeanConvFADC2PheVar < 0.)
    637     return -1.;
    638   return TMath::Sqrt(fMeanConvFADC2PheVar);
    639 }
     533
    640534
    641535// --------------------------------------------------------------------------
     
    803697  // Calculate the Error of Nphe
    804698  //
    805   const Float_t pheRelVar = ffactorsquareRelVar + meanSquareRelVar + rsigmaSquareRelVar;
    806   fPheFFactorMethodVar    = pheRelVar * fPheFFactorMethod * fPheFFactorMethod;
     699  const Float_t pheRelVar  = meanSquareRelVar + rsigmaSquareRelVar;
     700  fPheFFactorMethodStatVar = pheRelVar * fPheFFactorMethod * fPheFFactorMethod;
     701  fPheFFactorMethodSystVar = ffactorsquareRelVar * fPheFFactorMethod * fPheFFactorMethod;
    807702
    808703  if (IsDebug())
     
    816711    }
    817712
    818   if (fPheFFactorMethodVar < 0. )
     713  if (fPheFFactorMethodStatVar < 0. )
    819714    return kFALSE;
    820715
     
    879774  // the errors, but have to take account of this cancellation:
    880775  //
    881   Float_t convrelvar = ffactorsquareRelVar + GetMeanRelVar() + rsigmaSquareRelVar;
     776  Float_t convrelvar = GetMeanRelVar() + rsigmaSquareRelVar;
     777  if (IsHiGainSaturation())
     778      convrelvar += GetConversionHiLoRelVar();
     779
    882780  const Float_t limit = IsHiGainSaturation() ? fConvFFactorRelVarLimit * 4. : fConvFFactorRelVarLimit;
    883781
     
    906804    }
    907805 
    908   fMeanConvFADC2PheVar =  convrelvar * fMeanConvFADC2Phe * fMeanConvFADC2Phe;
     806  fMeanConvFADC2PheStatVar = convrelvar * fMeanConvFADC2Phe  * fMeanConvFADC2Phe;
     807  fMeanConvFADC2PheSystVar = ffactorsquareRelVar * fMeanConvFADC2Phe * fMeanConvFADC2Phe;
    909808 
    910809  SetFFactorMethodValid(kTRUE);
     
    967866                              + 0.25 * nphotonsrelvar;
    968867 
    969   fMeanFFactorFADC2PhotVar    = ffactorrelvar * fMeanFFactorFADC2Phot * fMeanFFactorFADC2Phot;
     868  fMeanFFactorFADC2PhotVar = ffactorrelvar * fMeanFFactorFADC2Phot * fMeanFFactorFADC2Phot;
    970869
    971870  if (IsDebug())
  • trunk/MagicSoft/Mars/mcalib/MCalibrationChargePix.h

    r7099 r7188  
    1717  static const Float_t fgConversionHiLo;         //! Default fConversionHiLo          (now set to: 10.)
    1818  static const Float_t fgConversionHiLoErr;      //! Default fConversionHiLoVar       (now set to: 2.5)
     19  static const Float_t fgConversionHiLoSigma;    //! Default fConversionHiLoSigma     (now set to: 2.5)
    1920  static const Float_t fgPheFFactorMethodLimit;  //! Default fPheFFactorMethodLimit   (now set to: 5.)
    2021  static const Float_t fgConvFFactorRelErrLimit; //! Default fConvFFactorRelErrLimit  (now set to: 0.35) 
     
    2526  Float_t fConversionHiLo;                  // Conversion factor betw. Hi Gain and Lo Gain 
    2627  Float_t fConversionHiLoVar;               // Variance Conversion factor betw. Hi and Lo Gain
     28  Float_t fConversionHiLoSigma;             // Sigma of conversion factor betw. Hi and Lo Gain
    2729  Float_t fConvFFactorRelVarLimit;          // Limit for acceptance rel. variance Conversion FADC2Phe
    2830  Float_t fLoGainPedRmsSquare;              // Pedestal RMS square of Low Gain
    2931  Float_t fLoGainPedRmsSquareVar;           // Pedestal RMS square Variance of Low Gain
    3032  Float_t fMeanConvFADC2Phe;                // Conversion factor (F-factor method)
    31   Float_t fMeanConvFADC2PheVar;             // Variance conversion factor (F-factor method)
     33  Float_t fMeanConvFADC2PheStatVar;         // Variance conversion factor, only stat. error
     34  Float_t fMeanConvFADC2PheSystVar;         // Variance conversion factor, only syst. error
    3235  Float_t fMeanFFactorFADC2Phot;            // Total mean F-Factor to photons (F-factor method)
    33   Float_t fMeanFFactorFADC2PhotVar;         // Variance mean F-Factor photons (F-factor method) 
     36  Float_t fMeanFFactorFADC2PhotVar;         // Variance mean F-Factor photons, only stat. error
    3437  Float_t fPed;                             // Pedestal (from MPedestalPix) times number FADC slices
    3538  Float_t fPedVar;                          // Variance of pedestal
    3639  Float_t fPedRms;                          // Pedestal RMS (from MPedestalPix) times sqrt nr. FADC slices
    3740  Float_t fPedRmsVar;                       // Pedestal RMS (from MPedestalPix) times sqrt nr. FADC slices
    38   Float_t fPheFFactorMethod;                // Number Phe's calculated (F-factor method)
    39   Float_t fPheFFactorMethodVar;             // Variance number of Phe's (F-factor method)
     41  Float_t fPheFFactorMethod;                // Number Phe's calculated  with F-factor method)
     42  Float_t fPheFFactorMethodStatVar;         // Variance number of Phe's, only stat. error
     43  Float_t fPheFFactorMethodSystVar;         // Variance number of Phe's, only syst. error
    4044  Float_t fPheFFactorMethodLimit;           // Min. number Photo-electrons for pix to be accepted.
    4145  Float_t fRSigmaSquare;                    // Square of Reduced sigma
     
    5155 
    5256public:
     57  MCalibrationChargePix(const char *name=NULL, const char *title=NULL);
    5358
    54   MCalibrationChargePix(const char *name=NULL, const char *title=NULL);
    55  
     59  // TObject
    5660  void Clear(Option_t *o="");
     61  void Copy(TObject& object) const;
    5762
    5863  // Setter
    59   void SetAbsTimeMean ( const Float_t f ) { fAbsTimeMean = f; }
    60   void SetAbsTimeRms  ( const Float_t f ) { fAbsTimeRms  = f; }
    61   void SetConversionHiLo    ( const Float_t c=fgConversionHiLo    )        { fConversionHiLo    = c;       }
    62   void SetConversionHiLoErr ( const Float_t e=fgConversionHiLoErr )        { fConversionHiLoVar = e*e;     }
     64  /*
    6365  void SetConvFFactorRelErrLimit   ( const Float_t f=fgConvFFactorRelErrLimit) { fConvFFactorRelVarLimit = f*f;}
    64   void SetFFactorMethodValid   ( const Bool_t b = kTRUE );
    6566  void SetMeanConvFADC2Phe      ( const Float_t f)                          { fMeanConvFADC2Phe       = f; }
    6667  void SetMeanConvFADC2PheVar   ( const Float_t f)                          { fMeanConvFADC2PheVar    = f; }
    6768  void SetMeanFFactorFADC2Phot  ( const Float_t f)                          { fMeanFFactorFADC2Phot   = f; }
    68   void SetPedestal              ( const Float_t ped, const Float_t pedrms, const Float_t pederr);
    69   void SetPed                   ( const Float_t ped, const Float_t pederr); 
    70   void SetPedRMS              ( const Float_t pedrms, const Float_t pedrmserr); 
    7169  void SetPheFFactorMethod      ( const Float_t f)                          { fPheFFactorMethod       = f; }
    7270  void SetPheFFactorMethodVar   ( const Float_t f)                          { fPheFFactorMethodVar    = f; } 
    7371  void SetPheFFactorMethodLimit ( const Float_t f=fgPheFFactorMethodLimit ) { fPheFFactorMethodLimit  = f; }
    7472  void SetNumSaturated          ( const Int_t   i)                          { fNumSaturated           = i; }
     73  */
     74  void SetFFactorMethodValid     (const Bool_t b = kTRUE );
     75  void SetPedestal               (const Float_t ped, const Float_t pedrms, const Float_t pederr);
     76  void SetPed                   ( const Float_t ped, const Float_t pederr); 
     77  void SetPedRMS              ( const Float_t pedrms, const Float_t pedrmserr); 
     78
     79  void SetAbsTimeMean            (const Float_t f)                          { fAbsTimeMean            = f; }
     80  void SetAbsTimeRms             (const Float_t f)                          { fAbsTimeRms             = f; }
     81  void SetConversionHiLo         (const Float_t c=fgConversionHiLo        ) { fConversionHiLo         = c; }
     82  void SetConversionHiLoErr      (const Float_t e=fgConversionHiLoErr     ) { fConversionHiLoVar      = e*e;}
     83  void SetConversionHiLoSigma    (const Float_t s=fgConversionHiLoSigma   ) { fConversionHiLoSigma    = s; }
     84  void SetConvFFactorRelErrLimit (const Float_t f=fgConvFFactorRelErrLimit) { fConvFFactorRelVarLimit = f*f;}
     85  void SetMeanConvFADC2Phe       (const Float_t f)                          { fMeanConvFADC2Phe       = f; }
     86  void SetMeanConvFADC2PheVar    (const Float_t f)                          { fMeanConvFADC2PheStatVar= f; }
     87  void SetMeanConvFADC2PheSystVar(const Float_t f)                          { fMeanConvFADC2PheSystVar= f; }
     88  void SetMeanFFactorFADC2Phot   (const Float_t f)                          { fMeanFFactorFADC2Phot   = f; }
     89  void SetNumSaturated           (const Int_t   i)                          { fNumSaturated           = i; }
     90  void SetPheFFactorMethod       (const Float_t f)                          { fPheFFactorMethod       = f; }
     91  void SetPheFFactorMethodVar    (const Float_t f)                          { fPheFFactorMethodStatVar= f; }
     92  void SetPheFFactorMethodSystVar(const Float_t f)                          { fPheFFactorMethodSystVar= f; }
     93  void SetPheFFactorMethodLimit  (const Float_t f=fgPheFFactorMethodLimit ) { fPheFFactorMethodLimit  = f; }
    7594 
    7695  // Getters
    77   Float_t GetAbsTimeMean             () const { return fAbsTimeMean;             }
    78   Float_t GetAbsTimeRms              () const { return fAbsTimeRms;              }
    79   Float_t GetConversionHiLo          () const { return fConversionHiLo;          }
    80   Float_t GetConversionHiLoErr       () const;
    81   Float_t GetConvertedMean           () const;
    82   Float_t GetConvertedMeanErr        () const;
    83   Float_t GetConvertedSigma          () const;
    84   Float_t GetConvertedSigmaErr       () const;
    85   Float_t GetConvertedRSigma         () const;
    86   Float_t GetConvertedRSigmaErr      () const;
    87   Float_t GetConvertedRSigmaSquare   () const; 
    88   Float_t GetMeanConvFADC2Phe        () const { return fMeanConvFADC2Phe;        }
    89   Float_t GetMeanConvFADC2PheErr     () const;
    90   Float_t GetMeanConvFADC2PheVar     () const { return fMeanConvFADC2PheVar;     }
    91   Float_t GetMeanFFactorFADC2Phot    () const { return fMeanFFactorFADC2Phot;    }
    92   Float_t GetMeanFFactorFADC2PhotErr () const;
    93   Float_t GetMeanFFactorFADC2PhotVar () const { return fMeanFFactorFADC2PhotVar; }   
    94   Int_t   GetNumSaturated            () const { return fNumSaturated;            }
    95   Float_t GetPed                     () const { return fPed;                     }
    96   Float_t GetPedErr                  () const { return TMath::Sqrt(fPedVar);     }
    97   Float_t GetPedRms                  () const;
    98   Float_t GetPedRmsErr               () const;
    99   Float_t GetPheFFactorMethod        () const { return fPheFFactorMethod;        }   
    100   Float_t GetPheFFactorMethodErr     () const;
    101   Float_t GetPheFFactorMethodVar     () const { return fPheFFactorMethodVar;     }
    102   Float_t GetPheFFactorMethodRelVar  () const;
    103   Float_t GetRSigma                  () const;
    104   Float_t GetRSigmaErr               () const;
    105   Float_t GetRSigmaRelVar            () const;
    106   Float_t GetRSigmaPerCharge         () const;
    107   Float_t GetRSigmaPerChargeErr      () const;
     96  Float_t GetAbsTimeMean        () const { return fAbsTimeMean;    }
     97  Float_t GetAbsTimeRms         () const { return fAbsTimeRms;     }
     98  Float_t GetConversionHiLo     () const { return fConversionHiLo; }
     99  Float_t GetConversionHiLoErr  () const { return fConversionHiLoVar<0 ? -1 : TMath::Sqrt(fConversionHiLoVar);  }
     100  Float_t GetConversionHiLoSigma() const { return fConversionHiLoSigma; }
     101  Float_t GetConvertedMean      () const { return IsHiGainSaturation() ? fLoGainMean * fConversionHiLo : fHiGainMean; }
     102  Float_t GetConvertedMeanErr   () const
     103  {
     104      if (!IsHiGainSaturation())
     105          return GetMeanErr();
     106      const Float_t logainrelvar = GetLoGainMeanRelVar();
     107      return logainrelvar<0 ? -1 : TMath::Sqrt(logainrelvar + GetConversionHiLoRelVar()) * GetConvertedMean();
     108  }
     109  Float_t GetConvertedSigma() const { return IsHiGainSaturation() ? fLoGainSigma * fConversionHiLo : fHiGainSigma; }
     110  Float_t GetConvertedSigmaErr() const
     111  {
     112    if (!IsHiGainSaturation())
     113        return GetSigmaErr();
    108114
    109   Bool_t IsFFactorMethodValid        () const;
     115    if (fLoGainSigmaVar<0 || fLoGainSigma<0)
     116        return -1.;
     117
     118    const Float_t sigmaRelVar = fLoGainSigmaVar/(fLoGainSigma*fLoGainSigma);
     119    return TMath::Sqrt(sigmaRelVar+GetConversionHiLoRelVar()) * GetConvertedSigma();
     120  }
     121  Float_t GetConvertedRSigma() const
     122  {
     123      if (fRSigmaSquare < 0)
     124          return -1;
     125      const Float_t rsigma = TMath::Sqrt(fRSigmaSquare);
     126      return IsHiGainSaturation() ? rsigma*fConversionHiLo : rsigma ;
     127  }
     128  Float_t GetConvertedRSigmaErr() const;
     129  Float_t GetConvertedRSigmaSquare() const
     130  {
     131      if (fRSigmaSquare < 0)
     132          return -1;
     133      return IsHiGainSaturation() ? fRSigmaSquare*fConversionHiLo*fConversionHiLo : fRSigmaSquare ;
     134  }
     135  Float_t GetMeanConvFADC2Phe() const { return fMeanConvFADC2Phe; }
     136  Float_t GetMeanConvFADC2PheErr  () const { return fMeanConvFADC2PheStatVar<0 ? -1 : TMath::Sqrt(fMeanConvFADC2PheStatVar); }
     137  Float_t GetMeanConvFADC2PheSystErr() const { return fMeanConvFADC2PheSystVar<0 ? -1 : TMath::Sqrt(fMeanConvFADC2PheSystVar); }
     138  Float_t GetMeanConvFADC2PheTotErr() const
     139  {
     140      if (fMeanConvFADC2PheSystVar<0 || fMeanConvFADC2PheStatVar<0)
     141          return -1.;
     142      return TMath::Sqrt(fMeanConvFADC2PheSystVar+fMeanConvFADC2PheStatVar);
     143  }
     144  Float_t GetFFactorFADC2Phe        () const { return gkFFactor;   }
     145  Float_t GetMeanConvFADC2PheVar    () const { return fMeanConvFADC2PheStatVar; }
     146  Float_t GetMeanConvFADC2PheSystVar() const { return fMeanConvFADC2PheSystVar; }
     147
     148  Float_t GetMeanFFactorFADC2Phot   () const { return fMeanFFactorFADC2Phot;    }
     149  Float_t GetMeanFFactorFADC2PhotErr() const { return fMeanFFactorFADC2PhotVar<0 ? -1. : TMath::Sqrt(fMeanFFactorFADC2PhotVar); }
     150  Float_t GetMeanFFactorFADC2PhotVar() const { return fMeanFFactorFADC2PhotVar; }
     151  Int_t   GetNumSaturated           () const { return fNumSaturated;            }   
     152  Float_t GetPed                    () const { return fPed;                     }
     153  Float_t GetPedErr                 () const { return TMath::Sqrt(fPedVar);     }
     154  Float_t GetPedRms                 () const
     155  {
     156    if (!IsHiGainSaturation())
     157        return fPedRms;
     158    return fLoGainPedRmsSquare<0 ? -1 : TMath::Sqrt(fLoGainPedRmsSquare);
     159  }
     160  Float_t GetPedRmsErr              () const;
     161  Float_t GetPheFFactorMethod       () const { return fPheFFactorMethod; }
     162  Float_t GetPheFFactorMethodErr    () const { return fPheFFactorMethodStatVar<0 ? -1 : TMath::Sqrt(fPheFFactorMethodStatVar); }
     163  Float_t GetPheFFactorMethodSystErr() const { return fPheFFactorMethodSystVar<0 ? -1 : TMath::Sqrt(fPheFFactorMethodSystVar); }
     164  Float_t GetPheFFactorMethodTotErr () const
     165  {
     166      if (fPheFFactorMethodStatVar<0 || fPheFFactorMethodSystVar<0)
     167          return -1.;
     168      return TMath::Sqrt(fPheFFactorMethodStatVar+fPheFFactorMethodSystVar);
     169  }
     170  Float_t GetPheFFactorMethodVar    () const { return fPheFFactorMethodStatVar; }
     171  Float_t GetPheFFactorMethodSystVar() const { return fPheFFactorMethodSystVar; }
     172  Float_t GetPheFFactorMethodRelVar () const { return fPheFFactorMethodStatVar<=0 ? -1 : fPheFFactorMethodStatVar / (fPheFFactorMethod * fPheFFactorMethod); }
     173  Float_t GetPheFFactorMethodRelSystVar() const { return fPheFFactorMethodSystVar<=0 ? -1. : fPheFFactorMethodSystVar / (fPheFFactorMethod * fPheFFactorMethod); }
     174  Float_t GetRSigma                 () const { return fRSigmaSquare<0 ? -1 : TMath::Sqrt(fRSigmaSquare); }
     175
     176  Float_t GetRSigmaErr         () const;
     177  Float_t GetRSigmaRelVar      () const;
     178  Float_t GetRSigmaPerCharge   () const;
     179  Float_t GetRSigmaPerChargeErr() const;
     180
     181  Bool_t IsFFactorMethodValid() const;
    110182
    111183  // Calculations
  • trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoPix.cc

    r5749 r7188  
    1818!   Author(s): Markus Gaug   02/2004 <mailto:markus@ifae.es>
    1919!
    20 !   Copyright: MAGIC Software Development, 2000-2004
     20!   Copyright: MAGIC Software Development, 2000-2005
    2121!
    2222\* ======================================================================== */
    2323/////////////////////////////////////////////////////////////////////////////
    24 //                                                                         //
    25 // MCalibrationHiLoPix                                                     //
    26 //                                                                         //
    27 // Storage container for high-gain vs. low-gain charge calibration results //
    28 // of one Pixel (PMT).                                                     //
    29 // The following "calibration" constants can be retrieved:                 //
     24//
     25// MCalibrationHiLoPix
     26//
     27// Storage container for high-gain vs. low-gain charge calibration results
     28// of one Pixel (PMT).
     29// The following "calibration" constants can be retrieved:
    3030// - GetHiLoRatio(): The mean conversion High-gain vs. Low-gain
    3131//   with which the low-gain result has to be multiplied
    32 // - GetHiLoSigma(): The Gauss sigma of histogrammed High-gain vs. Low-gain ratios
     32// - GetHiLoSigma(): The Gauss sigma of histogrammed High-gain vs.
     33//   Low-gain ratios
    3334//
    34 // See also: MHCalibrationHiLoPix, MHCalibrationHiLoCam              //
    35 //                                                                         //
     35// See also: MHCalibrationHiLoPix, MHCalibrationHiLoCam
     36//
     37// Class Version 2:
     38//  + Float_t fOffsetPerSlice;  // Offset from fit (per FADC slice)
     39//  + Float_t fGainRatio;       // Ratio of gains from fit
     40//
    3641/////////////////////////////////////////////////////////////////////////////
    3742#include "MCalibrationHiLoPix.h"
     
    4651//
    4752MCalibrationHiLoPix::MCalibrationHiLoPix(const char *name, const char *title)
     53  : fOffsetPerSlice(-9999.), fGainRatio(-1.)
    4854{
    4955
     
    5258
    5359}
    54 
  • trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoPix.h

    r5946 r7188  
    99{
    1010private:
     11    Float_t fOffsetPerSlice;            // Offset from fit (per FADC slice)
     12    Float_t fGainRatio;                 // Ratio of gains from fit
    1113
    1214public:
     15  MCalibrationHiLoPix(const char *name=NULL, const char *title=NULL);
    1316
    14   MCalibrationHiLoPix(const char *name=NULL, const char *title=NULL);
    15   ~MCalibrationHiLoPix() {}
    16  
     17  // Setter
     18  void SetGainRatio     (const Float_t f) { fGainRatio      = f; }
     19  void SetOffsetPerSlice(const Float_t f) { fOffsetPerSlice = f; }
     20
     21  // Getter
    1722  Float_t GetHiLoChargeRatio()         const { return GetHiGainMean();     }
    1823  Float_t GetHiLoChargeRatioErr()      const { return GetHiGainMeanErr();  }
     
    2631  Float_t GetHiLoTimeDiffSigmaErr()    const { return GetLoGainSigmaErr(); }
    2732  Float_t GetHiLoTimeDiffProb()        const { return GetLoGainProb();     }
     33  Float_t GetGainRatio     ()          const { return fGainRatio;          }
     34  Float_t GetOffsetPerSlice()          const { return fOffsetPerSlice;     }
    2835
    29   ClassDef(MCalibrationHiLoPix, 1)      // Container HiLo conversion Calibration Results Pixel
     36  ClassDef(MCalibrationHiLoPix, 2)      // Container HiLo conversion Calibration Results Pixel
    3037};
    3138
  • trunk/MagicSoft/Mars/mcalib/MCalibrationIntensityChargeCam.cc

    r6412 r7188  
    594594          varerr[i] = pix.GetRSigmaErr();
    595595        }
    596       if (option.Contains("abstime"))
     596      if (option.Contains("abstimemean"))
    597597        {
    598598          var   [i] = pix.GetAbsTimeMean();
    599599          varerr[i] = pix.GetAbsTimeRms();
    600600        }
     601      if (option.Contains("abstimerms"))
     602        {
     603          var   [i] = pix.GetAbsTimeRms();
     604          varerr[i] = pix.GetAbsTimeRms()/2.;
     605        }
    601606      if (option.Contains("blackout"))
    602607        {
     
    663668          var   [i] = pix.GetRSigmaPerCharge();
    664669          varerr[i] = pix.GetRSigmaPerChargeErr();
     670        }
     671      if (option.Contains("conversionfactor"))
     672        {
     673          const MCalibrationChargePix &apix = (MCalibrationChargePix&)cam->GetAverageArea(0);
     674          const Float_t mean = pix.GetConvertedMean();
     675          const Float_t phe  = apix.GetPheFFactorMethod();
     676
     677          var[i]    = phe/mean;
     678          varerr[i] = TMath::Sqrt(apix.GetPheFFactorMethodErr()*apix.GetPheFFactorMethodErr()/mean/mean
     679                                  + phe*phe/mean/mean/mean/mean*pix.GetConvertedMeanErr()*pix.GetConvertedMeanErr());
    665680        }
    666681  }
     
    735750         
    736751          if (option.Contains("rsigma"))
     752              pvar = pix.GetRSigma();
     753          if (option.Contains("abstimemean"))
     754              pvar = pix.GetAbsTimeMean();
     755          if (option.Contains("abstimerms"))
     756              pvar = pix.GetAbsTimeRms();
     757          if (option.Contains("conversionhilo"))
     758              pvar = pix.GetConversionHiLo();
     759          if (option.Contains("convertedmean"))
     760              pvar = pix.GetConvertedMean();
     761          if (option.Contains("convertedsigma"))
     762              pvar = pix.GetConvertedSigma();
     763          if (option.Contains("convertedrsigma"))
     764              pvar = pix.GetConvertedRSigma();
     765          if (option.Contains("meanconvfadc2phe"))
     766              pvar = pix.GetMeanConvFADC2Phe();
     767          if (option.Contains("meanffactorfadc2phot"))
     768              pvar = pix.GetMeanFFactorFADC2Phot();
     769          if (option.Contains("ped"))
     770              pvar = pix.GetPed();
     771          if (option.Contains("pedrms"))
     772              pvar = pix.GetPedRms();
     773          if (option.Contains("pheffactormethod"))
     774              pvar = pix.GetPheFFactorMethod();
     775          if (option.Contains("rsigmapercharge"))
     776              pvar = pix.GetRSigmaPerCharge();
     777          if (option.Contains("conversionfactor"))
     778          {
     779              const MCalibrationChargePix &apix = (MCalibrationChargePix&)cam->GetAverageArea(aidx);
     780              pvar = apix.GetPheFFactorMethod()/pix.GetConvertedMean();
     781          }
     782
     783
     784          variab  += pvar;
     785          variab2 += pvar*pvar;
     786          num++;
     787         
     788          camcharge.Fill(j,pvar);
     789          camcharge.SetUsed(j);
     790        }
     791     
     792      if (num > 1)
     793        {
     794          variab  /= num;
     795          variance = (variab2 - variab*variab*num) / (num-1);
     796
     797          vararea[i]    = variab;
     798          varareaerr[i] = variance>0 ? TMath::Sqrt(variance/num) : 999999999.;
     799
     800          //
     801          // Make also a Gauss-fit to the distributions. The RMS can be determined by
     802          // outlier, thus we look at the sigma and the RMS and take the smaller one, afterwards.
     803          //
     804          h = camcharge.ProjectionS(TArrayI(),TArrayI(1,&aidx),"_py",750);
     805          h->SetDirectory(NULL);
     806          h->Fit("gaus","QL");
     807          TF1 *fit = h->GetFunction("gaus");
     808
     809          Float_t ci2   = fit->GetChisquare();
     810          Float_t sigma = fit->GetParameter(2);
     811
     812          if (ci2 > 500. || sigma > varareaerr[i])
     813            {
     814              h->Fit("gaus","QLM");
     815              fit = h->GetFunction("gaus");
     816
     817              ci2   = fit->GetChisquare();
     818              sigma = fit->GetParameter(2);
     819            }
     820         
     821          const Float_t mean  = fit->GetParameter(1);
     822          const Float_t ndf   = fit->GetNDF();
     823         
     824          *fLog << inf << "Camera Nr: " << i << endl;
     825          *fLog << inf << option.Data() << " area idx: " << aidx << " Results: " << endl;
     826          *fLog << inf << "Mean: " << Form("%4.3f",mean)
     827                << "+-" << Form("%4.3f",fit->GetParError(1))
     828                << "  Sigma: " << Form("%4.3f",sigma) << "+-" << Form("%4.3f",fit->GetParError(2))
     829                << "  Chisquare: " << Form("%4.3f",ci2) << "  NDF  : " << ndf << endl;
     830          delete h;
     831          gROOT->GetListOfFunctions()->Remove(fit);
     832
     833          if (sigma<varareaerr[i] && ndf>2 && ci2<500.)
     834            {
     835              vararea   [i] = mean;
     836              varareaerr[i] = sigma/TMath::Sqrt((Float_t)num);
     837            }
     838        }
     839      else
     840        {
     841          vararea[i]    = -1.;
     842          varareaerr[i] = 0.;
     843        }
     844
     845      nr[i] = i;
     846      nrerr[i] = 0.;
     847    }
     848 
     849  TGraphErrors *gr = new TGraphErrors(size,
     850                                     nr.GetArray(),vararea.GetArray(),
     851                                     nrerr.GetArray(),varareaerr.GetArray());
     852  gr->SetTitle(Form("%s Area %3i Average",option.Data(),aidx));
     853  gr->GetXaxis()->SetTitle("Camera Nr.");
     854  //  gr->GetYaxis()->SetTitle("<Q> [1]");     
     855  return gr;
     856}
     857
     858
     859// -------------------------------------------------------------------
     860//
     861// Returns a TGraphErrors with the mean effective number of photon
     862// vs. the calibration camera number. With the string 'method', different
     863// calibration methods can be called.
     864//
     865TGraphErrors *MCalibrationIntensityChargeCam::GetPhotVsTime( const Option_t *method )
     866{
     867 
     868  const Int_t size = GetSize();
     869 
     870  if (size == 0)
     871    return NULL;
     872
     873  TString option(method);
     874
     875  TArrayF photarr(size);
     876  TArrayF photarrerr(size);
     877  TArrayF nr(size);
     878  TArrayF nrerr(size);
     879 
     880  for (Int_t i=0;i<GetSize();i++)
     881    {
     882      //
     883      // Get the calibration cam from the intensity cam
     884      //
     885      MCalibrationChargeCam *cam = (MCalibrationChargeCam*)GetCam(i);
     886
     887      //
     888      // Get the calibration pix from the calibration cam
     889      //
     890      Float_t phot    = 0.;
     891      Float_t photerr = 0.;
     892
     893      if (option.Contains("BlindPixel"))
     894        {
     895          phot    = cam->GetNumPhotonsBlindPixelMethod();
     896          photerr = cam->GetNumPhotonsBlindPixelMethodErr();
     897        }
     898      if (option.Contains("FFactor"))
     899        {
     900          phot    = cam->GetNumPhotonsFFactorMethod();
     901          photerr = cam->GetNumPhotonsFFactorMethodErr();
     902        }
     903      if (option.Contains("PINDiode"))
     904        {
     905          phot    = cam->GetNumPhotonsPINDiodeMethod();
     906          photerr = cam->GetNumPhotonsPINDiodeMethodErr();
     907        }
     908
     909      photarr[i]       = phot;
     910      photarrerr[i]    = photerr;
     911
     912      nr[i] = i;
     913      nrerr[i] = 0.;
     914    }
     915 
     916  TGraphErrors *gr = new TGraphErrors(size,
     917                                     nr.GetArray(),photarr.GetArray(),
     918                                     nrerr.GetArray(),photarrerr.GetArray());
     919  gr->SetTitle("Photons Average");
     920  gr->GetXaxis()->SetTitle("Camera Nr.");
     921  gr->GetYaxis()->SetTitle("<N_{phot}> [1]");
     922  return gr;
     923}
     924
     925// -------------------------------------------------------------------
     926//
     927// Returns a TGraphErrors with the mean effective number of photo-electrons per
     928// area index 'aidx' vs. the calibration camera number
     929//
     930TGraphErrors *MCalibrationIntensityChargeCam::GetPhePerAreaVsTime( const Int_t aidx, const MGeomCam &geom)
     931{
     932 
     933  const Int_t size = GetSize();
     934 
     935  if (size == 0)
     936    return NULL;
     937 
     938  TArrayF phearea(size);
     939  TArrayF pheareaerr(size);
     940  TArrayF time(size);
     941  TArrayF timeerr(size);
     942 
     943  for (Int_t i=0;i<GetSize();i++)
     944    {
     945      //
     946      // Get the calibration cam from the intensity cam
     947      //
     948      MCalibrationChargeCam *cam = (MCalibrationChargeCam*)GetCam(i);
     949
     950      //
     951      // Get the calibration pix from the calibration cam
     952      //
     953      const MCalibrationChargePix &apix = (MCalibrationChargePix&)cam->GetAverageArea(aidx);
     954      const Float_t phe          = apix.GetPheFFactorMethod();
     955      const Float_t pheerr       = apix.GetPheFFactorMethodErr();
     956
     957      phearea[i]       = phe;
     958      pheareaerr[i]    = pheerr;
     959
     960      time[i] = i;
     961      timeerr[i] = 0.;
     962    }
     963 
     964  TGraphErrors *gr = new TGraphErrors(size,
     965                                     time.GetArray(),phearea.GetArray(),
     966                                     timeerr.GetArray(),pheareaerr.GetArray());
     967  gr->SetTitle(Form("Phes Area %d Average",aidx));
     968  gr->GetXaxis()->SetTitle("Camera Nr.");
     969  gr->GetYaxis()->SetTitle("<N_{phe}> [1]");
     970  return gr;
     971}
     972
     973// -------------------------------------------------------------------
     974//
     975// Returns a TGraphErrors with the event-by-event averaged charge per
     976// area index 'aidx' vs. the calibration camera number
     977//
     978TGraphErrors *MCalibrationIntensityChargeCam::GetChargePerAreaVsTime( const Int_t aidx, const MGeomCam &geom)
     979{
     980 
     981  const Int_t size = GetSize();
     982 
     983  if (size == 0)
     984    return NULL;
     985 
     986  TArrayF chargearea(size);
     987  TArrayF chargeareaerr(size);
     988  TArrayF nr(size);
     989  TArrayF nrerr(size);
     990 
     991  for (Int_t i=0;i<GetSize();i++)
     992    {
     993      //
     994      // Get the calibration cam from the intensity cam
     995      //
     996      MCalibrationChargeCam *cam = (MCalibrationChargeCam*)GetCam(i);
     997
     998      //
     999      // Get the calibration pix from the calibration cam
     1000      //
     1001      const MCalibrationChargePix &apix = (MCalibrationChargePix&)cam->GetAverageArea(aidx);
     1002      const Float_t charge          = apix.GetConvertedMean();
     1003      const Float_t chargeerr       = apix.GetConvertedSigma();
     1004
     1005      chargearea[i]       = charge;
     1006      chargeareaerr[i]    = chargeerr;
     1007
     1008      nr[i]    = i;
     1009      nrerr[i] = 0.;
     1010    }
     1011 
     1012  TGraphErrors *gr = new TGraphErrors(size,
     1013                                     nr.GetArray(),chargearea.GetArray(),
     1014                                     nrerr.GetArray(),chargeareaerr.GetArray());
     1015  gr->SetTitle(Form("Averaged Charges Area Idx %d",aidx));
     1016  gr->GetXaxis()->SetTitle("Camera Nr.");
     1017  gr->GetYaxis()->SetTitle("<Q> [FADC cnts]");     
     1018  return gr;
     1019}
     1020
     1021TH1F *MCalibrationIntensityChargeCam::GetVarFluctuations( const Int_t aidx, const MGeomCam &geom, const Option_t *varname )
     1022{
     1023 
     1024  const Int_t size = GetSize();
     1025 
     1026  if (size == 0)
     1027    return NULL;
     1028 
     1029  TString option(varname);
     1030  option.ToLower();
     1031 
     1032  TH1F *hist = new TH1F("hist",Form("%s - Rel. Fluctuations %s Pixel",option.Data(),aidx ? "Outer" : "Inner"),
     1033                        200,0.,100.);
     1034  hist->SetXTitle("Relative Fluctuation [%]");
     1035  hist->SetYTitle("Nr. channels [1]"); 
     1036  hist->SetFillColor(kRed+aidx);
     1037
     1038  MCalibrationChargeCam *cam = (MCalibrationChargeCam*)GetCam();
     1039
     1040  //
     1041  // Loop over pixels
     1042  //
     1043  for (Int_t npix=0;npix<cam->GetSize();npix++)
     1044    {
     1045      if (geom[npix].GetAidx() != aidx)
     1046        continue;
     1047
     1048      Double_t variab   = 0.;
     1049      Double_t variab2  = 0.;
     1050      Double_t variance = 0.;
     1051      Int_t    num      = 0;
     1052      Float_t  pvar     = 0.;
     1053      Float_t  relrms   = 99.9;
     1054      //
     1055      // Loop over the Cams for each pixel
     1056      //
     1057      for (Int_t i=0; i<GetSize(); i++)
     1058        {
     1059          MCalibrationChargeCam *cam = (MCalibrationChargeCam*)GetCam(i);
     1060          //
     1061          // Get the calibration pix from the calibration cam
     1062          //
     1063          MCalibrationChargePix &pix = (MCalibrationChargePix&)(*cam)[npix];
     1064          //
     1065          // Don't use bad pixels
     1066          //
     1067          if (!pix.IsFFactorMethodValid())
     1068            continue;
     1069
     1070          if (option.Contains("rsigma"))
    7371071            pvar = pix.GetRSigma();
    738           if (option.Contains("abstime"))
     1072          if (option.Contains("abstimemean"))
    7391073            pvar = pix.GetAbsTimeMean();
     1074          if (option.Contains("abstimerms"))
     1075            pvar = pix.GetAbsTimeRms();
    7401076          if (option.Contains("conversionhilo"))
    7411077            pvar = pix.GetConversionHiLo();
     
    7581094          if (option.Contains("rsigmapercharge"))
    7591095            pvar = pix.GetRSigmaPerCharge();
     1096          if (option.Contains("conversionfactor"))
     1097          {
     1098              const MCalibrationChargePix &apix = (MCalibrationChargePix&)cam->GetAverageArea(0);
     1099              pvar = apix.GetPheFFactorMethod()/pix.GetConvertedMean();
     1100          }
     1101
    7601102
    7611103          variab  += pvar;
    7621104          variab2 += pvar*pvar;
    7631105          num++;
    764          
    765           camcharge.Fill(j,pvar);
    766           camcharge.SetUsed(j);
    767         }
    768      
    769       if (num > 1)
    770         {
    771           variab  /= num;
    772           variance = (variab2 - variab*variab*num) / (num-1);
    773 
    774           vararea[i] = variab;
    775           if (variance > 0.)
    776             varareaerr[i] = TMath::Sqrt(variance);
    777           else
    778             varareaerr[i] = 999999999.;
    779 
    780           //
    781           // Make also a Gauss-fit to the distributions. The RMS can be determined by
    782           // outlier, thus we look at the sigma and the RMS and take the smaller one, afterwards.
    783           //
    784           h = camcharge.ProjectionS(TArrayI(),TArrayI(1,&aidx),"_py",750);
    785           h->SetDirectory(NULL);
    786           h->Fit("gaus","QL");
    787           TF1 *fit = h->GetFunction("gaus");
    788 
    789           Float_t ci2   = fit->GetChisquare();
    790           Float_t sigma = fit->GetParameter(2);
    791 
    792           if (ci2 > 500. || sigma > varareaerr[i])
    793             {
    794               h->Fit("gaus","QLM");
    795               fit = h->GetFunction("gaus");
    796 
    797               ci2   = fit->GetChisquare();
    798               sigma = fit->GetParameter(2);
    799             }
    800          
    801           const Float_t mean  = fit->GetParameter(1);
    802           const Float_t ndf   = fit->GetNDF();
    803          
    804           *fLog << inf << "Camera Nr: " << i << endl;
    805           *fLog << inf << option.Data() << " area idx: " << aidx << " Results: " << endl;
    806           *fLog << inf << "Mean: " << Form("%4.3f",mean)
    807                 << "+-" << Form("%4.3f",fit->GetParError(1))
    808                 << "  Sigma: " << Form("%4.3f",sigma) << "+-" << Form("%4.3f",fit->GetParError(2))
    809                 << "  Chisquare: " << Form("%4.3f",fit->GetChisquare()) << "  NDF  : " << ndf << endl;         
    810           delete h;
    811           gROOT->GetListOfFunctions()->Remove(fit);
    812 
    813           if (sigma < varareaerr[i] && ndf > 2)
    814             {
    815               vararea   [i] = mean;
    816               varareaerr[i] = sigma;
    817             }
    818         }
    819       else
    820         {
    821           vararea[i]    = -1.;
    822           varareaerr[i] = 0.;
    823         }
    824 
    825       nr[i] = i;
    826       nrerr[i] = 0.;
    827     }
    828  
    829   TGraphErrors *gr = new TGraphErrors(size,
    830                                      nr.GetArray(),vararea.GetArray(),
    831                                      nrerr.GetArray(),varareaerr.GetArray());
    832   gr->SetTitle(Form("%s Area %3i Average",option.Data(),aidx));
    833   gr->GetXaxis()->SetTitle("Camera Nr.");
    834   //  gr->GetYaxis()->SetTitle("<Q> [1]");     
    835   return gr;
    836 }
    837 
    838 
    839 // -------------------------------------------------------------------
    840 //
    841 // Returns a TGraphErrors with the mean effective number of photon
    842 // vs. the calibration camera number. With the string 'method', different
    843 // calibration methods can be called.
    844 //
    845 TGraphErrors *MCalibrationIntensityChargeCam::GetPhotVsTime( const Option_t *method )
    846 {
    847  
    848   const Int_t size = GetSize();
    849  
    850   if (size == 0)
    851     return NULL;
    852 
    853   TString option(method);
    854 
    855   TArrayF photarr(size);
    856   TArrayF photarrerr(size);
    857   TArrayF nr(size);
    858   TArrayF nrerr(size);
    859  
    860   for (Int_t i=0;i<GetSize();i++)
    861     {
    862       //
    863       // Get the calibration cam from the intensity cam
    864       //
    865       MCalibrationChargeCam *cam = (MCalibrationChargeCam*)GetCam(i);
    866 
    867       //
    868       // Get the calibration pix from the calibration cam
    869       //
    870       Float_t phot    = 0.;
    871       Float_t photerr = 0.;
    872 
    873       if (option.Contains("BlindPixel"))
    874         {
    875           phot    = cam->GetNumPhotonsBlindPixelMethod();
    876           photerr = cam->GetNumPhotonsBlindPixelMethodErr();
    877         }
    878       if (option.Contains("FFactor"))
    879         {
    880           phot    = cam->GetNumPhotonsFFactorMethod();
    881           photerr = cam->GetNumPhotonsFFactorMethodErr();
    882         }
    883       if (option.Contains("PINDiode"))
    884         {
    885           phot    = cam->GetNumPhotonsPINDiodeMethod();
    886           photerr = cam->GetNumPhotonsPINDiodeMethodErr();
    887         }
    888 
    889       photarr[i]       = phot;
    890       photarrerr[i]    = photerr;
    891 
    892       nr[i] = i;
    893       nrerr[i] = 0.;
    894     }
    895  
    896   TGraphErrors *gr = new TGraphErrors(size,
    897                                      nr.GetArray(),photarr.GetArray(),
    898                                      nrerr.GetArray(),photarrerr.GetArray());
    899   gr->SetTitle("Photons Average");
    900   gr->GetXaxis()->SetTitle("Camera Nr.");
    901   gr->GetYaxis()->SetTitle("<N_phot> [1]");     
    902   return gr;
    903 }
    904 
    905 // -------------------------------------------------------------------
    906 //
    907 // Returns a TGraphErrors with the mean effective number of photo-electrons per
    908 // area index 'aidx' vs. the calibration camera number
    909 //
    910 TGraphErrors *MCalibrationIntensityChargeCam::GetPhePerAreaVsTime( const Int_t aidx, const MGeomCam &geom)
    911 {
    912  
    913   const Int_t size = GetSize();
    914  
    915   if (size == 0)
    916     return NULL;
    917  
    918   TArrayF phearea(size);
    919   TArrayF pheareaerr(size);
    920   TArrayF time(size);
    921   TArrayF timeerr(size);
    922  
    923   for (Int_t i=0;i<GetSize();i++)
    924     {
    925       //
    926       // Get the calibration cam from the intensity cam
    927       //
    928       MCalibrationChargeCam *cam = (MCalibrationChargeCam*)GetCam(i);
    929 
    930       //
    931       // Get the calibration pix from the calibration cam
    932       //
    933       const MCalibrationChargePix &apix = (MCalibrationChargePix&)cam->GetAverageArea(aidx);
    934       const Float_t phe          = apix.GetPheFFactorMethod();
    935       const Float_t pheerr       = apix.GetPheFFactorMethodErr();
    936 
    937       phearea[i]       = phe;
    938       pheareaerr[i]    = pheerr;
    939 
    940       time[i] = i;
    941       timeerr[i] = 0.;
    942     }
    943  
    944   TGraphErrors *gr = new TGraphErrors(size,
    945                                      time.GetArray(),phearea.GetArray(),
    946                                      timeerr.GetArray(),pheareaerr.GetArray());
    947   gr->SetTitle(Form("Phes Area %d Average",aidx));
    948   gr->GetXaxis()->SetTitle("Camera Nr.");
    949   gr->GetYaxis()->SetTitle("<N_phes> [1]");     
    950   return gr;
    951 }
    952 
    953 // -------------------------------------------------------------------
    954 //
    955 // Returns a TGraphErrors with the event-by-event averaged charge per
    956 // area index 'aidx' vs. the calibration camera number
    957 //
    958 TGraphErrors *MCalibrationIntensityChargeCam::GetChargePerAreaVsTime( const Int_t aidx, const MGeomCam &geom)
    959 {
    960  
    961   const Int_t size = GetSize();
    962  
    963   if (size == 0)
    964     return NULL;
    965  
    966   TArrayF chargearea(size);
    967   TArrayF chargeareaerr(size);
    968   TArrayF nr(size);
    969   TArrayF nrerr(size);
    970  
    971   for (Int_t i=0;i<GetSize();i++)
    972     {
    973       //
    974       // Get the calibration cam from the intensity cam
    975       //
    976       MCalibrationChargeCam *cam = (MCalibrationChargeCam*)GetCam(i);
    977 
    978       //
    979       // Get the calibration pix from the calibration cam
    980       //
    981       const MCalibrationChargePix &apix = (MCalibrationChargePix&)cam->GetAverageArea(aidx);
    982       const Float_t charge          = apix.GetConvertedMean();
    983       const Float_t chargeerr       = apix.GetConvertedSigma();
    984 
    985       chargearea[i]       = charge;
    986       chargeareaerr[i]    = chargeerr;
    987 
    988       nr[i]    = i;
    989       nrerr[i] = 0.;
    990     }
    991  
    992   TGraphErrors *gr = new TGraphErrors(size,
    993                                      nr.GetArray(),chargearea.GetArray(),
    994                                      nrerr.GetArray(),chargeareaerr.GetArray());
    995   gr->SetTitle(Form("Averaged Charges Area Idx %d",aidx));
    996   gr->GetXaxis()->SetTitle("Camera Nr.");
    997   gr->GetYaxis()->SetTitle("<Q> [FADC cnts]");     
    998   return gr;
    999 }
    1000 
    1001 TH1F *MCalibrationIntensityChargeCam::GetVarFluctuations( const Int_t aidx, const MGeomCam &geom, const Option_t *varname )
    1002 {
    1003  
    1004   const Int_t size = GetSize();
    1005  
    1006   if (size == 0)
    1007     return NULL;
    1008  
    1009   TString option(varname);
    1010  
    1011   TH1F *hist = new TH1F("hist",Form("%s - Rel. Fluctuations %s Pixel",option.Data(),aidx ? "Outer" : "Inner"),
    1012                         200,0.,100.);
    1013   hist->SetXTitle("Relative Fluctuation [%]");
    1014   hist->SetYTitle("Nr. channels [1]"); 
    1015   hist->SetFillColor(kRed+aidx);
    1016 
    1017   MCalibrationChargeCam *cam = (MCalibrationChargeCam*)GetCam();
    1018 
    1019   //
    1020   // Loop over pixels
    1021   //
    1022   for (Int_t npix=0;npix<cam->GetSize();npix++)
    1023     {
    1024       if (geom[npix].GetAidx() != aidx)
    1025         continue;
    1026 
    1027       Double_t variab   = 0.;
    1028       Double_t variab2  = 0.;
    1029       Double_t variance = 0.;
    1030       Int_t    num      = 0;
    1031       Float_t  pvar     = 0.;
    1032       Float_t  relrms   = 99.9;
    1033       //
    1034       // Loop over the Cams for each pixel
    1035       //
    1036       for (Int_t i=0; i<GetSize(); i++)
    1037         {
    1038           MCalibrationChargeCam *cam = (MCalibrationChargeCam*)GetCam(i);
    1039           //
    1040           // Get the calibration pix from the calibration cam
    1041           //
    1042           MCalibrationChargePix &pix = (MCalibrationChargePix&)(*cam)[npix];
    1043           //
    1044           // Don't use bad pixels
    1045           //
    1046           if (!pix.IsFFactorMethodValid())
    1047             continue;
    1048 
    1049           if (option.Contains("RSigma"))
    1050             pvar = pix.GetRSigma();
    1051           if (option.Contains("AbsTime"))
    1052             pvar = pix.GetAbsTimeMean();
    1053           if (option.Contains("ConversionHiLo"))
    1054             pvar = pix.GetConversionHiLo();
    1055           if (option.Contains("ConvertedMean"))
    1056             pvar = pix.GetConvertedMean();
    1057           if (option.Contains("ConvertedSigma"))
    1058             pvar = pix.GetConvertedSigma();
    1059           if (option.Contains("ConvertedRSigma"))
    1060             pvar = pix.GetConvertedRSigma();
    1061           if (option.Contains("MeanConvFADC2Phe"))
    1062             pvar = pix.GetMeanConvFADC2Phe();
    1063           if (option.Contains("MeanFFactorFADC2Phot"))
    1064             pvar = pix.GetMeanFFactorFADC2Phot();
    1065           if (option.Contains("Ped"))
    1066             pvar = pix.GetPed();
    1067           if (option.Contains("PedRms"))
    1068             pvar = pix.GetPedRms();
    1069           if (option.Contains("PheFFactorMethod"))
    1070             pvar = pix.GetPheFFactorMethod();
    1071           if (option.Contains("RSigmaPerCharge"))
    1072             pvar = pix.GetRSigmaPerCharge();
    1073 
    1074           variab  += pvar;
    1075           variab2 += pvar*pvar;
    1076           num++;
    10771106        }
    10781107
  • trunk/MagicSoft/Mars/mcalib/MCalibrationRelTimeCalc.cc

    r6963 r7188  
    4040//                  - FinalizeBadPixels()
    4141//
     42// Class Version 2:
     43//  + Byte_t fCheckFlags; // Bit-field to hold the possible check flags
     44//
     45//
    4246//  Input Containers:
    4347//   MCalibrationRelTimeCam
     
    7175#include "MBadPixelsPix.h"
    7276
    73 
    7477ClassImp(MCalibrationRelTimeCalc);
    7578
     
    7780
    7881const Float_t MCalibrationRelTimeCalc::fgRelTimeResolutionLimit = 1.0;
     82
    7983// --------------------------------------------------------------------------
    8084//
     
    97101  fName  = name  ? name  : "MCalibrationRelTimeCalc";
    98102  fTitle = title ? title : "Task to finalize the relative time calibration";
    99  
     103
     104  SetCheckFitResults       ( kFALSE );
     105  SetCheckDeviatingBehavior( kFALSE );
     106  SetCheckHistOverflow     ( kFALSE );
     107  SetCheckOscillations     ( kFALSE );
     108
    100109  SetRelTimeResolutionLimit();
    101110  SetOutputPath();
    102111  SetOutputFile("");
    103  
     112
    104113  Clear();
    105  
     114
    106115}
    107116
     
    268277
    269278  PrintUncalibrated(MBadPixelsPix::kDeviatingTimeResolution,   
    270                     Form("%s%2.1f%s","Time resolution less than ",fRelTimeResolutionLimit," FADC slices from Mean:   "));
     279                    Form("%s%2.1f%s","Time resol. less than ",fRelTimeResolutionLimit," FADC sl. from Mean:     "));
    271280  PrintUncalibrated(MBadPixelsPix::kRelTimeOscillating,   
    272281                    "Pixels with changing Rel. Times over time:             ");
     
    437446      MCalibrationRelTimePix &pix = (MCalibrationRelTimePix&)(*relcam)[i];
    438447
    439       if (bad.IsUncalibrated( MBadPixelsPix::kDeviatingTimeResolution))
    440         bad.SetUnsuitable(   MBadPixelsPix::kUnreliableRun   );
    441  
    442       if (bad.IsUncalibrated( MBadPixelsPix::kRelTimeNotFitted))
    443         bad.SetUnsuitable(   MBadPixelsPix::kUnreliableRun   );
    444  
    445       if (bad.IsUncalibrated( MBadPixelsPix::kRelTimeOscillating))
    446         bad.SetUnsuitable(   MBadPixelsPix::kUnreliableRun   );
    447  
    448       if (bad.IsUnsuitable(   MBadPixelsPix::kUnsuitableRun    ))
    449         pix.SetExcluded();
    450 
    451     }
     448      if (IsCheckDeviatingBehavior())
     449          if (bad.IsUncalibrated(MBadPixelsPix::kDeviatingTimeResolution))
     450              bad.SetUnsuitable(MBadPixelsPix::kUnreliableRun);
     451
     452      if (IsCheckFitResults())
     453          if (bad.IsUncalibrated(MBadPixelsPix::kRelTimeNotFitted))
     454              bad.SetUnsuitable(MBadPixelsPix::kUnreliableRun);
     455
     456      if (IsCheckOscillations())
     457          if (bad.IsUncalibrated(MBadPixelsPix::kRelTimeOscillating))
     458              bad.SetUnsuitable(MBadPixelsPix::kUnreliableRun);
     459
     460      if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
     461          pix.SetExcluded();
     462    }
     463
    452464}
    453465
     
    488500
    489501  if (fGeom->InheritsFrom("MGeomCamMagic"))
    490     *fLog << " " << setw(7) << "Uncalibrated Pixels:            "
    491           << Form("%s%3i%s%3i","Inner: ",counts[0]," Outer: ",counts[1]) << endl;
     502    *fLog << " " << setw(7) << "Uncalibrated Pixels:            Inner: "
     503          << Form("%3i",counts[0]) << " Outer: " << Form("%3i",counts[1]) << endl;
    492504
    493505  counts.Reset();
     
    506518    relcam->SetNumUnreliable(counts[aidx], aidx);
    507519
    508   *fLog << " " << setw(7) << "Unreliable Pixels:              "
    509         << Form("%s%3i%s%3i","Inner: ",counts[0]," Outer: ",counts[1]) << endl;
     520  *fLog << " " << setw(7) << "Unreliable Pixels:              Inner: "
     521        << Form("%3i",counts[0]) << " Outer: " << Form("%3i",counts[1]) << endl;
    510522
    511523}
     
    545557void MCalibrationRelTimeCalc::SetOutputPath(TString path)
    546558{
    547   fOutputPath = path;
    548   if (fOutputPath.EndsWith("/"))
    549     fOutputPath = fOutputPath(0, fOutputPath.Length()-1);
     559    fOutputPath = path;
     560    if (fOutputPath.EndsWith("/"))
     561        fOutputPath = fOutputPath(0, fOutputPath.Length()-1);
    550562}
    551563
     
    556568const char* MCalibrationRelTimeCalc::GetOutputFile()
    557569{
    558   return Form("%s/%s", (const char*)fOutputPath, (const char*)fOutputFile);
    559 }
     570    return Form("%s/%s", (const char*)fOutputPath, (const char*)fOutputFile);
     571}
     572
     573
     574// --------------------------------------------------------------------------
     575//
     576// MCalibrationRelTimeCam.CheckFitResults: Yes
     577// MCalibrationRelTimeCam.CheckDeviatingBehavior: Yes
     578// MCalibrationRelTimeCam.CheckHistOverflow: Yes
     579// MCalibrationRelTimeCam.CheckOscillations: Yes
     580//
     581Int_t MCalibrationRelTimeCalc::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
     582{
     583    Bool_t rc = kFALSE;
     584
     585    if (IsEnvDefined(env, prefix, "CheckFitResults", print))
     586    {
     587        SetCheckFitResults(GetEnvValue(env, prefix, "CheckFitResults", IsCheckFitResults()));
     588        rc = kTRUE;
     589    }
     590    if (IsEnvDefined(env, prefix, "CheckDeviatingBehavior", print))
     591    {
     592        SetCheckDeviatingBehavior(GetEnvValue(env, prefix, "CheckDeviatingBehavior", IsCheckDeviatingBehavior()));
     593        rc = kTRUE;
     594    }
     595    if (IsEnvDefined(env, prefix, "CheckHistOverflow", print))
     596    {
     597        SetCheckHistOverflow(GetEnvValue(env, prefix, "CheckHistOverflow", IsCheckHistOverflow()));
     598        rc = kTRUE;
     599    }
     600
     601    if (IsEnvDefined(env, prefix, "CheckOscillations", print))
     602    {
     603        SetCheckOscillations(GetEnvValue(env, prefix, "CheckOscillations", IsCheckOscillations()));
     604        rc = kTRUE;
     605    }
     606
     607    return rc;
     608}
  • trunk/MagicSoft/Mars/mcalib/MCalibrationRelTimeCalc.h

    r5679 r7188  
    4343  MGeomCam                   *fGeom;             //! Camera geometry
    4444
     45  // enums
     46  enum  Check_t
     47  {
     48      kCheckHistOverflow,
     49      kCheckDeviatingBehavior,
     50      kCheckOscillations,
     51      kCheckFitResults
     52  };                                             // Possible Checks
     53
     54  Byte_t fCheckFlags;                            // Bit-field to hold the possible check flags
     55
    4556  enum  { kDebug };                              //  Possible flags
    4657
     
    5667  void   PrintUncalibrated( MBadPixelsPix::UncalibratedType_t typ, const char *text) const;
    5768
     69  // Query checks
     70  Bool_t IsCheckDeviatingBehavior() const { return TESTBIT(fCheckFlags,kCheckDeviatingBehavior); }
     71  Bool_t IsCheckHistOverflow     () const { return TESTBIT(fCheckFlags,kCheckHistOverflow);      }
     72  Bool_t IsCheckOscillations     () const { return TESTBIT(fCheckFlags,kCheckOscillations);      }
     73  Bool_t IsCheckFitResults       () const { return TESTBIT(fCheckFlags,kCheckFitResults);        }
     74
     75  // MTask
    5876  Bool_t ReInit     (MParList *pList);
    5977  Int_t  Process    () { return kTRUE; }
    6078  Int_t  PostProcess();
    6179
     80  // MParContainer
     81  Int_t  ReadEnv(const TEnv &env, TString prefix, Bool_t print);
     82
    6283public:
    63 
    6484  MCalibrationRelTimeCalc(const char *name=NULL, const char *title=NULL);
    6585
     86  // TObject
    6687  void Clear(const Option_t *o="");
    67  
     88
     89  // MCalibrationRelTimeCalc
    6890  Int_t Finalize();
    6991 
     92  // Getter
    7093  Bool_t IsDebug() const  {  return TESTBIT(fFlags,kDebug); }
    7194
    72   void SetDebug                 ( const Bool_t  b=kTRUE ) { b ? SETBIT(fFlags,kDebug) : CLRBIT(fFlags,kDebug); }   
     95  // Setter
    7396  void SetOutputPath            ( TString path="."                         );
    7497  void SetOutputFile            ( TString file="TimeCalibStat.txt"         ) { fOutputFile          = file; }
    7598  void SetRelTimeResolutionLimit( const Float_t f=fgRelTimeResolutionLimit ) { fRelTimeResolutionLimit = f; }
    7699
    77   ClassDef(MCalibrationRelTimeCalc, 1)   // Task finalizing the relative time Calibration
     100  // Checks
     101  void SetCheckFitResults(const Bool_t b=kTRUE)        { b ? SETBIT(fCheckFlags,kCheckFitResults) : CLRBIT(fCheckFlags,kCheckFitResults); }
     102  void SetCheckDeviatingBehavior(const Bool_t b=kTRUE) { b ? SETBIT(fCheckFlags,kCheckDeviatingBehavior) : CLRBIT(fCheckFlags,kCheckDeviatingBehavior); }
     103  void SetCheckHistOverflow(const Bool_t b=kTRUE)      { b ? SETBIT(fCheckFlags,kCheckHistOverflow) : CLRBIT(fCheckFlags,kCheckHistOverflow); }
     104  void SetCheckOscillations(const Bool_t b=kTRUE)      { b ? SETBIT(fCheckFlags,kCheckOscillations) : CLRBIT(fCheckFlags,kCheckOscillations); }
     105  void SetDebug(const Bool_t b=kTRUE)                  { b ? SETBIT(fFlags, kDebug) : CLRBIT(fFlags, kDebug); }
     106
     107  ClassDef(MCalibrationRelTimeCalc, 2)   // Task finalizing the relative time Calibration
    78108};
    79109
  • trunk/MagicSoft/Mars/mcalib/MMcCalibrationCalc.cc

    r7005 r7188  
    7171using namespace std;
    7272
    73 MMcCalibrationCalc::MMcCalibrationCalc(const char *name, const char *title)
     73// --------------------------------------------------------------------------
     74//
     75// Constructor. Default value for fMinSize is 1000 ADC counts. This must be
     76// set in general by the user (SetMinSize), since it depends among other things
     77// on the signal extractor used.
     78//
     79MMcCalibrationCalc::MMcCalibrationCalc(const char *name, const char *title): fMinSize(1000)
    7480{
    7581    fName  = name  ? name  : "MMcCalibrationCalc";
    7682    fTitle = title ? title : "Calculate and write conversion factors into MCalibrationChargeCam and MCalibrationQECam containers";
    77 
    78     fHistADC2PhotEl = new TH1F(AddSerialNumber("ADC2PhotEl"), "log10(fPhotElfromShower/fSize)", 1500, -3., 3.);
    79     fHistADC2PhotEl->SetXTitle("log_{10}(fPhotElfromShower / fSize) [photel/ADC count]");
    80 
    81 
    82     fHistPhot2PhotEl = new TH1F(AddSerialNumber("Phot2PhotEl"), "Photon conversion efficiency", 1000, 0., 1.);
    83     fHistPhot2PhotEl->SetXTitle("Overall photon conversion efficiency [photoelectron/photon]");
    84 
    8583}
    8684
     
    9391Bool_t MMcCalibrationCalc::CheckRunType(MParList *pList) const
    9492{
    95     const MRawRunHeader *run = (MRawRunHeader*)pList->FindObject("MRawRunHeader");
     93    const MRawRunHeader *run = (MRawRunHeader*)pList->FindObject(AddSerialNumber("MRawRunHeader"));
    9694    if (!run)
    9795    {
    98         *fLog << warn << "Warning - cannot check file type, MRawRunHeader not found." << endl;
    99         return kTRUE;
     96      *fLog << warn << "Warning - cannot check file type, " << AddSerialNumber("MRawRunHeader")
     97            << " not found." << endl;
     98      return kTRUE;
    10099    }
    101100
     
    105104// --------------------------------------------------------------------------
    106105//
    107 // Make sure, that there is an MCalibrationCam Object in the Parameter List.
     106// Look for all necessary containers and create histograms
    108107//
    109108Int_t MMcCalibrationCalc::PreProcess(MParList *pList)
    110109{
    111     fHistADC2PhotEl->Reset();
    112     fHistPhot2PhotEl->Reset();
    113 
    114110    fADC2PhotEl = 0;
    115111    fPhot2PhotEl = 0;
     
    156152        return kFALSE;
    157153    }
     154
     155    //
     156    // Create histograms:
     157    //
     158
     159    fHistADC2PhotEl = new TH1F(AddSerialNumber("ADC2PhotEl"), "log10(fPhotElfromShower/fSize)", 1500, -3., 3.);
     160    fHistADC2PhotEl->SetXTitle("log_{10}(fPhotElfromShower / fSize) [photel/ADC count]");
     161
     162
     163    fHistPhot2PhotEl = new TH1F(AddSerialNumber("Phot2PhotEl"), "Photon conversion efficiency", 1000, 0., 1.);
     164    fHistPhot2PhotEl->SetXTitle("Overall photon conversion efficiency [photoelectron/photon]");
     165
    158166
    159167    return kTRUE;
     
    210218  // perpendicular to the camera plane.
    211219  //
    212   // FIXME! We should look for AddSerialNumber("MMcConfigRunHeader") but
    213   // for the moment the stereo version of camera does not write one such
    214   // header per telescope (it should!)
    215   //
    216   MMcConfigRunHeader* mcconfig = (MMcConfigRunHeader*) pList->FindObject("MMcConfigRunHeader");
     220  MMcConfigRunHeader* mcconfig = (MMcConfigRunHeader*) pList->FindObject(AddSerialNumber("MMcConfigRunHeader"));
    217221  if (!mcconfig)
    218222    {
    219       *fLog << err << "MMcConfigRunHeader" << " not found... aborting." << endl;
     223      *fLog << err << AddSerialNumber("MMcConfigRunHeader") << " not found... aborting." << endl;
    220224      return kFALSE;
    221225    }
     
    253257    //
    254258    // Exclude events with low Size (larger fluctuations)
    255     // FIXME? The present cut (1000 "inner-pixel-counts") is somehow
    256     // arbitrary. Might it be optimized?
    257259    //   
    258260
    259     if (size < 1000)
     261    if (size < fMinSize)
    260262        return kTRUE;
    261263
     
    288290    }
    289291
    290     fPhot2PhotEl = fHistPhot2PhotEl->GetMean();   // Average quantum efficiency
     292    fPhot2PhotEl = fHistPhot2PhotEl->GetMean();   
     293    // Average quantum efficiency. For now we will set this value for all pixels, although
     294    // even in MC there may be differences from pixel to pixel, if the .dat file containing
     295    // QE vs lambda, input of the camera simulation, has different QE curves for different
     296    // pixels. FIXME?
     297
     298    MCalibrationQEPix &avqepix = (MCalibrationQEPix&)(fQECam->GetAverageArea(0));
     299    avqepix.SetAverageQE(fPhot2PhotEl); // Needed by MCalibrateData!
    291300
    292301    //
     
    335344        // FIXME: we are now assuming that all inner pixels have the same gain, and all
    336345        // outer pixels have the same gain (different from inner ones though). This can
    337         // only be like this in camera 0.7, but may change in future versions of camera.
     346        // only be like this in camera 0.7, but might change in future versions of camera.
    338347        //
    339348
     
    347356    }
    348357
     358
    349359    return kTRUE;
    350360}
  • trunk/MagicSoft/Mars/mcalib/MMcCalibrationCalc.h

    r4710 r7188  
    3535                                         // outer pixels w.r.t inner ones
    3636
     37    Float_t fMinSize;
     38    // Minimum SIZE (before calibration, ADC counts) an event must have to be considered in the
     39    // calculation of the calibration constants.
     40
    3741    TH1F*   fHistADC2PhotEl;
    3842    TH1F*   fHistPhot2PhotEl; // Histograms for monitoring the calibration.
     
    5054    TH1F*   GetHistPhot2PhotEl() { return fHistPhot2PhotEl; }
    5155
     56    void SetMinSize(Float_t x) { fMinSize = x; }
     57
    5258    ClassDef(MMcCalibrationCalc, 0)   // Task which obtains, for MC files, the calibration factor from ADC counts to photons.
    5359};
  • trunk/MagicSoft/Mars/mcalib/Makefile

    r6662 r7188  
    4040           MCalibrateRelTimes.cc \
    4141           MCalibrationIntensityCam.cc \
     42           MCalibrationIntensityConstCam.cc \
    4243           MCalibrationIntensityChargeCam.cc \
    4344           MCalibrationIntensityBlindCam.cc \
Note: See TracChangeset for help on using the changeset viewer.