Changeset 3650 for trunk


Ignore:
Timestamp:
04/05/04 16:44:56 (21 years ago)
Author:
gaug
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars
Files:
1 added
5 edited

Legend:

Unmodified
Added
Removed

  • trunk/MagicSoft/Mars/Changelog

    r3645 r3650  
    1919                                                 -*-*- END OF LINE -*-*-
    2020
     21 2004/04/05: Markus Gaug
     22
     23   * mcalib/MCalibrationPedCam.[h,cc]
     24   * mcalib/Makefile
     25   * mcalib/CalibLinkDef.h
     26     - new class to store information obtained from MHPedestalCam
     27
     28
    2129 2004/04/02: Markus Gaug
    2230
     31   * mcalib/MCalibrationChargePix.[h,cc]
    2332   * mcalib/MCalibrationPix.[h,cc]
    2433   * mcalib/MCalibrationCam.[h,cc]

  • trunk/MagicSoft/Mars/mcalib/MCalibrationChargePix.cc

    r3645 r3650  
    2222\* ======================================================================== */
    2323/////////////////////////////////////////////////////////////////////////////
    24 //                                                                         //
    25 // MCalibrationChargePix                                                   //
    26 //                                                                         //
    27 // Storage container to hold informations about the calibration values     //
    28 // values of one Pixel (PMT).                                              //
    29 //                                                                         //
     24//                                                                         
     25// MCalibrationChargePix                                                   
     26//                                                                         
     27// Storage container of the calibrated Charge of one pixel.
     28//                                                                         
    3029// The following values are initialized to meaningful values:
    3130//
    … …  
    3534//   with the Munich definition of the F-Factor, thus:
    3635//   F = Sigma(Out)/Mean(Out) * Mean(In)/Sigma(In)
    37 //   Mean F-Factor  = 1.15
    38 //   Error F-Factor = 0.02
    39 //
    40 // - Average QE: (email David Paneque, 14.2.04):
    41 //
    42 //  The conversion factor that comes purely from QE folded to a Cherenkov
    43 //  spectrum has to be multiplied by:
    44 //  * Plexiglass window -->> 0.96 X 0.96
    45 //  * PMT photoelectron collection efficiency -->> 0.9
    46 //  * Light guides efficiency -->> 0.94
    47 //
    48 //  Concerning the light guides effiency estimation... Daniel Ferenc
    49 //  is preparing some work (simulations) to estimate it. Yet so far, he has
    50 //  been busy with other stuff, and this work is still UNfinished.
    51 //
    52 //  The estimation I did comes from:
    53 //  1) Reflectivity of light guide walls is 85 % (aluminum)
    54 //  2) At ZERO degree light incidence, 37% of the light hits such walls
    55 //    (0.15X37%= 5.6% of light lost)
    56 //  3) When increasing the light incidence angle, more and more light hits
    57 //     the walls.
    58 //
    59 //  However, the loses due to larger amount of photons hitting the walls is more
    60 //  or less counteracted by the fact that more and more photon trajectories cross
    61 //  the PMT photocathode twice, increasing the effective sensitivity of the PMT.
    62 //
    63 //  Jurgen Gebauer did some quick measurements about this issue. I attach a
    64 //  plot. You can see that the angular dependence is (more or less) in agreement
    65 //  with a CosTheta function (below 20-25 degrees),
    66 //  which is the variation of teh entrance window cross section. So, in
    67 //  first approximation, no loses when increasing light incidence angle;
    68 //  and therefore, the factor 0.94.
    69 //
    70 //  So, summarizing... I would propose the following conversion factors
    71 //  (while working with CT1 cal box) in order to get the final number of photons
    72 //  from the detected measured size in ADC counts.
    73 //
    74 //  Nph = ADC * FmethodConversionFactor / ConvPhe-PhFactor
    75 //
    76 //  FmethodConversionFactor ; measured for individual pmts
    77 //
    78 //  ConvPhe-PhFactor = 0.98 * 0.23 * 0.90 * 0.94 * 0.96 * 0.96 = 0.18
    79 //
    80 //  I would not apply any smearing of this factor (which we have in nature),
    81 //  since we might be applying it to PMTs in the totally wrong direction.
    82 //
     36//   Mean F-Factor (gkFFactor)     = 1.15
     37//   Error F-Factor (gkFFactorErr) = 0.02
     38//
     39// The following variables are calculated inside this class:
     40// -  fLoGainPedRmsSquare and LoGainPedRmsSquareVar (see CalcLoGainPedestal())
     41// -  fRSigmaSquare and fRSigmaSquareVar            (see CalcReducedSigma() )
     42// -  fPheFFactorMethod and fPheFFactorMethodVar    see CalcFFactorMethod() )
     43//
     44// The following variables are set by MHCalibrationChargeCam:
     45// -  fAbsTimeMean and fAbsTimeRms
     46// -  all variables in MCalibrationPix
     47//
     48// The following variables are set by MCalibrationChargeCalc:
     49// - fPed, fPedVar and fPedRms                         
     50// - fMeanConversionFFactorMethod, fMeanConversionBlindPixelMethod,
     51//   fMeanConversionPINDiodeMethod and fMeanConversionCombinedMethod   
     52// - fConversionFFactorMethodVar, fConversionBlindPixelMethodVar
     53//   fConversionPINDiodeMethodVar and fConversionCombinedMethodVar
     54// - fSigmaConversionFFactorMethodm, fSigmaConversionBlindPixelMethod
     55//   fSigmaConversionPINDiodeMethod and fSigmaConversionCombinedMethod 
     56// - fTotalFFactorFFactorMethod, fTotalFFactorBlindPixelMethod   
     57//   fTotalFFactorPINDiodeMethod and fTotalFFactorCombinedMethod     
     58// - fTotalFFactorFFactorMethodVar, fTotalFFactorBlindPixelMethodVar,
     59//   fTotalFFactorPINDiodeMethodVar and fTotalFFactorCombinedMethodVar 
     60//
     61// The following variables are not yet implemented:
     62// - fConversionHiLo and fConversionHiLoVar (now set fixed to 10. +- 2.5)
    8363//
    8464//  Error of all variables are calculated by error-propagation. Note that internally,
    8565//  all error variables contain Variances in order to save the CPU-intensive square rooting
    8666//
     67// See also: MCalibrationChargeCam, MCalibrationChargeCalc,
     68//           MHCalibrationChargeCam, MHCalibrationChargePix
     69//
    8770/////////////////////////////////////////////////////////////////////////////
    8871#include "MCalibrationChargePix.h"
    … …  
    10992// Default Constructor:
    11093//
     94// Sets:
     95// - fCalibFlags to 0
     96// - fConversionHiLo to fgConversionHiLo
     97// - fConversionHiLoErr to fgConversionHiLoErr
     98// - PheFFactorMethodLimit to fgPheFFactorMethodLimit
     99//
     100// Calls:
     101// - Clear()
     102//
    111103MCalibrationChargePix::MCalibrationChargePix(const char *name, const char *title)
    112104    : fCalibFlags(0)
    … …  
    131123// ------------------------------------------------------------------------
    132124//
    133 // Invalidate values
     125// Sets:
     126// - all flags to 0
     127// - all variables to -1.
     128//
     129// Calls:
     130// - MCalibrationPix::Clear()
    134131//
    135132void MCalibrationChargePix::Clear(Option_t *o)
    … …  
    141138  SetCombinedMethodValid    ( kFALSE );
    142139
    143   fRSigma                           =  -1.;
    144   fRSigmaVar                        =  -1.;
     140  fRSigmaSquare                     =  -1.;
     141  fRSigmaSquareVar                  =  -1.;
    145142 
    146143  fPed                              =  -1.;
    … …  
    148145  fPedVar                           =  -1.;
    149146
    150   fLoGainPedRms                     =  -1.;
    151   fLoGainPedRmsVar                  =  -1.;
     147  fLoGainPedRmsSquare               =  -1.;
     148  fLoGainPedRmsSquareVar            =  -1.;
    152149
    153150  fAbsTimeMean                      =  -1.;
    … …  
    183180// --------------------------------------------------------------------------
    184181//
    185 // Set the pedestals from outside
    186 //
    187 void MCalibrationChargePix::SetPedestal(const Float_t ped, const Float_t pedrms, const Float_t pederr)
    188 {
    189 
    190   fPed       = ped;   
    191   fPedRms    = pedrms;
    192   fPedVar    = pederr*pederr;
    193 }
    194  
    195 
    196 
    197 // --------------------------------------------------------------------------
    198 //
    199 // Set the conversion factors from outside (only for MC)
     182// Set the conversion factors Blind Pixel Method from outside (only for MC)
     183//
     184void MCalibrationChargePix::SetConversionBlindPixelMethod(Float_t c, Float_t err, Float_t sig)
     185{
     186  fMeanConversionBlindPixelMethod  = c;
     187  fConversionBlindPixelMethodVar   = err*err;
     188  fSigmaConversionBlindPixelMethod = sig;
     189}
     190
     191
     192// --------------------------------------------------------------------------
     193//
     194// Set the conversion factors Combined Method from outside (only for MC)
     195//
     196void MCalibrationChargePix::SetConversionCombinedMethod(Float_t c, Float_t err, Float_t sig)
     197{
     198  fMeanConversionCombinedMethod  = c;
     199  fConversionCombinedMethodVar   = err*err;
     200  fSigmaConversionCombinedMethod = sig;
     201}
     202
     203
     204// --------------------------------------------------------------------------
     205//
     206// Set the conversion factors F-Factor Method from outside (only for MC)
    200207//
    201208void MCalibrationChargePix::SetConversionFFactorMethod(Float_t c, Float_t err, Float_t sig)
    … …  
    208215// --------------------------------------------------------------------------
    209216//
    210 // Set the conversion factors from outside (only for MC)
    211 //
    212 void MCalibrationChargePix::SetConversionCombinedMethod(Float_t c, Float_t err, Float_t sig)
    213 {
    214   fMeanConversionCombinedMethod  = c;
    215   fConversionCombinedMethodVar   = err*err;
    216   fSigmaConversionCombinedMethod = sig;
    217 }
    218 
    219 
    220 // --------------------------------------------------------------------------
    221 //
    222 // Set the conversion factors from outside (only for MC)
    223 //
    224 void MCalibrationChargePix::SetConversionBlindPixelMethod(Float_t c, Float_t err, Float_t sig)
    225 {
    226   fMeanConversionBlindPixelMethod  = c;
    227   fConversionBlindPixelMethodVar   = err*err;
    228   fSigmaConversionBlindPixelMethod = sig;
    229 }
    230 
    231 // --------------------------------------------------------------------------
    232 //
    233 // Set the conversion factors from outside (only for MC)
     217// Set the conversion factors PIN Diode Method from outside (only for MC)
    234218//
    235219void MCalibrationChargePix::SetConversionPINDiodeMethod(Float_t c, Float_t err, Float_t sig)
    … …  
    242226// --------------------------------------------------------------------------
    243227//
    244 // Set the Excluded Bit from outside
     228// Set the BlindPixelMethod Validity Bit from outside
    245229//
    246230void MCalibrationChargePix::SetBlindPixelMethodValid(const Bool_t b )
    … …  
    251235// --------------------------------------------------------------------------
    252236//
    253 // Set the Excluded Bit from outside
     237// Set the CombinedMethod Validity Bit from outside
     238//
     239void MCalibrationChargePix::SetCombinedMethodValid(const Bool_t b ) 
     240{
     241  b ?  SETBIT(fCalibFlags, kCombinedMethodValid) : CLRBIT(fCalibFlags, kCombinedMethodValid);
     242}
     243
     244// --------------------------------------------------------------------------
     245//
     246// Set the FFactorMethod Validity Bit from outside
    254247//
    255248void MCalibrationChargePix::SetFFactorMethodValid(const Bool_t b )
    … …  
    260253// --------------------------------------------------------------------------
    261254//
    262 // Set the Excluded Bit from outside
     255// Set the PINDiodeMethod Validity Bit from outside
    263256//
    264257void MCalibrationChargePix::SetPINDiodeMethodValid(const Bool_t b ) 
    … …  
    269262// --------------------------------------------------------------------------
    270263//
    271 // Set the Excluded Bit from outside
    272 //
    273 void MCalibrationChargePix::SetCombinedMethodValid(const Bool_t b ) 
    274 {
    275   b ?  SETBIT(fCalibFlags, kCombinedMethodValid) : CLRBIT(fCalibFlags, kCombinedMethodValid);
    276 }
    277 
    278 
     264// Set the pedestals from outside
     265//
     266void MCalibrationChargePix::SetPedestal(const Float_t ped, const Float_t pedrms, const Float_t pederr)
     267{
     268
     269  fPed       = ped;   
     270  fPedRms    = pedrms;
     271  fPedVar    = pederr*pederr;
     272}
     273 
     274
     275// --------------------------------------------------------------------------
     276//
     277// Get the Error of the Mean pedestals:
     278// Returns square root of fPedVar
     279//
     280Float_t  MCalibrationChargePix::GetPedErr()  const
     281{
     282  return  TMath::Sqrt(fPedVar);
     283}
     284
     285// --------------------------------------------------------------------------
     286//
     287// Get the pedestals RMS:
     288// - Test bit kHiGainSaturation:
     289//   If yes, return square root of fLoGainPedRmsSquare (if greater than 0, otherwise -1.),
     290//   If no,  return fPedRms
     291//
    279292Float_t  MCalibrationChargePix::GetPedRms()  const
    280293{
    281   return IsHiGainSaturation() ? fLoGainPedRms : fPedRms;
    282 }
    283 
     294
     295  if (IsHiGainSaturation())
     296    if (fLoGainPedRmsSquare < 0.)
     297      return -1.;
     298    else
     299      return TMath::Sqrt(fLoGainPedRmsSquare);
     300 
     301  return fPedRms;
     302}
     303
     304// --------------------------------------------------------------------------
     305//
     306// Get the Error of the pedestals RMS:
     307// - Test bit kHiGainSaturation:
     308//   If yes, return kLoGainPedRms, else fPedRms
     309//
    284310Float_t  MCalibrationChargePix::GetPedRmsErr()  const
    285311{
    286   return IsHiGainSaturation() ? TMath::Sqrt(fLoGainPedRmsVar) : TMath::Sqrt(fPedVar)/2.;
    287 }
    288 
    289 Float_t  MCalibrationChargePix::GetPedErr()  const
    290 {
    291   return  TMath::Sqrt(fPedVar);
    292 }
    293 
    294 
     312  if (IsHiGainSaturation())
     313    if (fLoGainPedRmsSquareVar < 0.)
     314      return -1.;
     315    else
     316      return TMath::Sqrt(0.25*fLoGainPedRmsSquareVar/fLoGainPedRmsSquare);
     317  else
     318    if (fPedVar < 0.)
     319      return -1.;
     320    else
     321      return TMath::Sqrt(fPedVar)/2.;
     322}
     323
     324
     325// --------------------------------------------------------------------------
     326//
     327// Get the Low Gain Mean:
     328// Returns fLoGainMean multiplied with fConversionHiLo
     329//
    295330Float_t MCalibrationChargePix::GetLoGainMean()  const
    296331{
    … …  
    298333}
    299334
     335// --------------------------------------------------------------------------
     336//
     337// Get the Error of the Low Gain Mean:
     338// Returns the quadratic sum of the relative low Gain Mean error and the
     339// the relative conversion High-to-Low error, mulitplied with GetLoGainMean()
     340//
    300341Float_t MCalibrationChargePix::GetLoGainMeanErr()  const
    301342{
    … …  
    304345                                 /( fLoGainMean * fLoGainMean );
    305346
    306   const Float_t conversionRelVar =  fConversionHiLoVar
    307                                  /( fConversionHiLo   * fConversionHiLo   );
    308 
    309   return TMath::Sqrt(chargeRelVar+conversionRelVar) * GetLoGainMean();
    310 }
    311 
     347  return TMath::Sqrt(chargeRelVar+GetConversionHiLoRelVar()) * GetLoGainMean();
     348}
     349
     350// --------------------------------------------------------------------------
     351//
     352// Get the Low Gain Sigma:
     353// Returns fLoGainSigma multiplied with fConversionHiLo
     354//
    312355Float_t MCalibrationChargePix::GetLoGainSigma()  const
    313356{
    … …  
    315358}
    316359
     360// --------------------------------------------------------------------------
     361//
     362// Get the Error of the Low Gain Sigma:
     363// Returns the quadratic sum of the relative low Gain Sigma error and the
     364// the relative conversion High-to-Low error, mulitplied with GetLoGainSigma()
     365//
    317366Float_t MCalibrationChargePix::GetLoGainSigmaErr()  const
    318367{
    … …  
    321370                                /( fLoGainSigma * fLoGainSigma );
    322371
    323   const Float_t conversionRelVar =  fConversionHiLoVar
    324                                  /( fConversionHiLo   * fConversionHiLo   );
    325 
    326   return TMath::Sqrt(sigmaRelVar+conversionRelVar) * GetLoGainSigma();
    327 }
    328 
     372  return TMath::Sqrt(sigmaRelVar+GetConversionHiLoRelVar()) * GetLoGainSigma();
     373}
     374
     375// --------------------------------------------------------------------------
     376//
     377// Get the reduced Sigma:
     378// - If fRSigmaSquare is smaller than 0 (i.e. has not yet been set), return -1.
     379// - Test bit kHiGainSaturation:
     380//   If yes, return square root of fRSigmaSquare, multiplied with fConversionHiLo,
     381//   If no , return square root of fRSigmaSquare
     382//
    329383Float_t MCalibrationChargePix::GetRSigma()  const
    330384{
    331   return IsHiGainSaturation() ? fRSigma*fConversionHiLo : fRSigma ;
     385  if (fRSigmaSquare < 0)
     386    return -1;
     387
     388  const Float_t rsigma = TMath::Sqrt(fRSigmaSquare);
     389 
     390  return IsHiGainSaturation() ? rsigma*fConversionHiLo : rsigma ;
    332391}
    333392
     393// --------------------------------------------------------------------------
     394//
     395// Get the error of the reduced Sigma:
     396// - If fRSigmaSquareVar is smaller than 0 (i.e. has not yet been set), return -1.
     397// - Calculate the absolute variance of the reduced sigma with the formula:
     398//   sigma variance = 0.25 * fRSigmaSquareVar / fRSigmaSquare
     399// - Test bit kHiGainSaturation:
     400//   If yes, returns the quadratic sum of the relative reduced Sigma error and the
     401//    the relative conversion High-to-Low error, mulitplied with GetRSigma()
     402//   Else returns the square root of rel. (0.25*fRSigmaSquareVar / fRSigmaSquare)
     403//
    334404Float_t MCalibrationChargePix::GetRSigmaErr()  const
    335405{
    336  if (IsHiGainSaturation())
    337     {
    338       const Float_t rsigmaRelVar     =  fRSigmaVar
    339                                     /( fRSigma * fRSigma );
    340       const Float_t conversionRelVar =  fConversionHiLoVar
    341                                      /( fConversionHiLo   * fConversionHiLo   );
    342       return TMath::Sqrt(rsigmaRelVar+conversionRelVar) * GetRSigma();
    343     }
     406
     407  if (fRSigmaSquareVar < 0)
     408    return -1;
     409
     410  //
     411  // SigmaSquareVar = 4. * Sigma * Sigma * Var(sigma)
     412  // ==> Var(sigma) = 0.25 * SigmaSquareVar / (Sigma * Sigma)
     413  //
     414  const Float_t rsigmaVar = 0.25 * fRSigmaSquareVar / fRSigmaSquare;
     415
     416  if (IsHiGainSaturation())
     417    return TMath::Sqrt(rsigmaVar/fRSigmaSquare + GetConversionHiLoRelVar()) * GetRSigma();
    344418 else
    345    return TMath::Sqrt(fRSigmaVar);
    346 
    347 }
    348 
     419   return TMath::Sqrt(rsigmaVar);
     420
     421}
     422
     423// --------------------------------------------------------------------------
     424//
     425// Get the conversion Error Hi-Gain to Low-Gain:
     426// - If fConversionHiLoVar is smaller than 0 (i.e. has not yet been set), return -1.
     427// - Else returns the square root of fConversionHiLoVar
     428//
    349429Float_t MCalibrationChargePix::GetConversionHiLoErr()  const
    350430{
    351431  if (fConversionHiLoVar < 0.)
    352432    return -1.;
     433
    353434  return TMath::Sqrt(fConversionHiLoVar);
    354435}
    355436
     437// --------------------------------------------------------------------------
     438//
     439// Get the error on the number of photo-electrons (F-Factor Method):
     440// - If fPheFFactorMethodVar is smaller than 0 (i.e. has not yet been set), return -1.
     441// - Else returns the square root of fPheFFactorMethodVar
     442//
    356443Float_t MCalibrationChargePix::GetPheFFactorMethodErr()  const
    357444{
    … …  
    361448}
    362449
     450// --------------------------------------------------------------------------
     451//
     452// Get the error on the mean conversion factor (Combined Method):
     453// - If fConversionCombinedMethodVar is smaller than 0 (i.e. has not yet been set), return -1.
     454// - Else returns the square root of fConversionCombinedMethodVar
     455//
    363456Float_t MCalibrationChargePix::GetConversionCombinedMethodErr()  const
    364457{
    … …  
    368461}
    369462
     463// --------------------------------------------------------------------------
     464//
     465// Get the error on the mean conversion factor (PINDiode Method):
     466// - If fConversionPINDiodeMethodVar is smaller than 0 (i.e. has not yet been set), return -1.
     467// - Else returns the square root of fConversionPINDiodeMethodVar
     468//
    370469Float_t MCalibrationChargePix::GetConversionPINDiodeMethodErr()  const
    371470{
    … …  
    375474}
    376475
     476// --------------------------------------------------------------------------
     477//
     478// Get the error on the mean conversion factor (BlindPixel  Method):
     479// - If fConversionBlindPixelMethodVar is smaller than 0 (i.e. has not yet been set), return -1.
     480// - Else returns the square root of fConversionBlindPixelMethodVar
     481//
    377482Float_t MCalibrationChargePix::GetConversionBlindPixelMethodErr()  const
    378483{
    … …  
    382487}
    383488
     489// --------------------------------------------------------------------------
     490//
     491// Get the error on the mean conversion factor (FFactor  Method):
     492// - If fConversionFFactorMethodVar is smaller than 0 (i.e. has not yet been set), return -1.
     493// - Else returns the square root of fConversionFFactorMethodVar
     494//
    384495Float_t MCalibrationChargePix::GetConversionFFactorMethodErr()  const
    385496{
    … …  
    389500}
    390501
     502// --------------------------------------------------------------------------
     503//
     504// Get the error on the total F-Factor (Combined  Method):
     505// - If fTotalFFactorCombinedMethodVar is smaller than 0 (i.e. has not yet been set), return -1.
     506// - Else returns the square root of fTotalFFactorCombinedMethodVar
     507//
    391508Float_t MCalibrationChargePix::GetTotalFFactorCombinedMethodErr()  const
    392509{
    … …  
    396513}
    397514
     515// --------------------------------------------------------------------------
     516//
     517// Get the error on the total F-Factor (PIN Diode  Method):
     518// - If fTotalPINDiodeCombinedMethodVar is smaller than 0 (i.e. has not yet been set), return -1.
     519// - Else returns the square root of fTotalPINDiodeCombinedMethodVar
     520//
    398521Float_t MCalibrationChargePix::GetTotalFFactorPINDiodeMethodErr()  const
    399522{
    … …  
    403526}
    404527
     528// --------------------------------------------------------------------------
     529//
     530// Get the error on the total F-Factor (Blind Pixel  Method):
     531// - If fTotalBlindPixelCombinedMethodVar is smaller than 0 (i.e. has not yet been set), return -1.
     532// - Else returns the square root of fTotalBlindPixelCombinedMethodVar
     533//
    405534Float_t MCalibrationChargePix::GetTotalFFactorBlindPixelMethodErr()  const
    406535{
    … …  
    410539}
    411540
     541// --------------------------------------------------------------------------
     542//
     543// Get the error on the total F-Factor (F-Factor  Method):
     544// - If fTotalFFactorCombinedMethodVar is smaller than 0 (i.e. has not yet been set), return -1.
     545// - Else returns the square root of fTotalFFactorCombinedMethodVar
     546//
    412547Float_t MCalibrationChargePix::GetTotalFFactorFFactorMethodErr()  const
    413548{
    … …  
    417552}
    418553
     554// --------------------------------------------------------------------------
     555//
     556// Test bit kBlindPixelMethodValid
     557//
    419558Bool_t MCalibrationChargePix::IsBlindPixelMethodValid() const
    420559{
    … …  
    422561}
    423562
     563// --------------------------------------------------------------------------
     564//
     565// Test bit kFFactorMethodValid
     566//
    424567Bool_t MCalibrationChargePix::IsFFactorMethodValid()   const
    425568{
    … …  
    427570}
    428571
     572// --------------------------------------------------------------------------
     573//
     574// Test bit kPINDiodeMethodValid
     575//
    429576Bool_t MCalibrationChargePix::IsPINDiodeMethodValid()  const
    430577{
    … …  
    432579}
    433580
     581// --------------------------------------------------------------------------
     582//
     583// Test bit kCombinedMethodValid
     584//
    434585Bool_t MCalibrationChargePix::IsCombinedMethodValid()  const
    435586{
    … …  
    437588}
    438589
    439 //
     590// ----------------------------------------------------------------------------
    440591//
     592// - If fSigma  is smaller than 0 (i.e. has not yet been set), return kFALSE
     593// - If fPedRms is smaller than 0 (i.e. has not yet been set), return kFALSE
     594//
     595// Calculate the reduced sigma:
     596// - Test bit IsHiGainSaturation() for the Sigma:
     597//   If yes, take fLoGainSigma and fLoGainSigmaVar
     598//   If no , take fHiGainSigma and fHiGainSigmaVar
     599//
     600// - Test bit IsHiGainSaturation() for the pedRMS:
     601//   If yes, take fLoGainPedRmsSquare and fLoGainPedRmsSquareVar
     602//   If no , take fPedRms and fPedRmsVar
     603//
     604// - Calculate the reduced sigma with the formula:
     605//   fRSigmaSquare = Sigma*Sigma - pedRMS*pedRMS
     606//
     607// - Calculate the variance of the reduced sigma with the formula:
     608//   fRSigmaSquareVar = 4.* (sigmaVar*Sigma*Sigma + pedRmsVar*pedRMS*pedRMS)
    441609//
    442610Bool_t MCalibrationChargePix::CalcReducedSigma()
    443611{
     612
     613  if (GetSigma() < 0.)
     614    return kFALSE;
     615 
     616  if (GetPedRms() < 0.)
     617    return kFALSE;
    444618
    445619  const Float_t sigma    = IsHiGainSaturation() ? fLoGainSigma    : fHiGainSigma   ;
    446620  const Float_t sigmavar = IsHiGainSaturation() ? fLoGainSigmaVar : fHiGainSigmaVar;
    447 
    448   const Float_t sigmaSquare     =      sigma     * sigma;
    449   const Float_t sigmaSquareVar  = 4.*  sigmavar  * sigmaSquare;
    450 
    451   Float_t pedRmsSquare ;         
    452   Float_t pedRmsSquareVar;
    453  
    454   if (IsHiGainSaturation())
    455     {
    456       pedRmsSquare     =      fLoGainPedRms    * fLoGainPedRms;
    457       pedRmsSquareVar  =  4.* fLoGainPedRmsVar * pedRmsSquare;
    458     }
    459   else
    460   {
    461       pedRmsSquare       =  fPedRms * fPedRms;                                         
    462       pedRmsSquareVar    =  fPedVar * pedRmsSquare; // fPedRmsErr = fPedErr/2.
    463   }
     621  const Float_t pedRmsSquare    = IsHiGainSaturation() ? fLoGainPedRmsSquare    : fPedRms*fPedRms;
     622  const Float_t pedRmsSquareVar = IsHiGainSaturation() ? fLoGainPedRmsSquareVar : 0.25*fPedVar*pedRmsSquare;
     623
     624  const Float_t sigmaSquare    = sigma     * sigma;
     625  const Float_t sigmaSquareVar = sigmavar  * sigmaSquare;
     626
    464627  //
    465628  // Calculate the reduced sigmas
    466629  //
    467   const Float_t rsigmasquare = sigmaSquare - pedRmsSquare;
    468   if (rsigmasquare <= 0.)
     630  fRSigmaSquare = sigmaSquare - pedRmsSquare;
     631  if (fRSigmaSquare <= 0.)
    469632    {
    470633      *fLog << warn
    … …  
    473636      return kFALSE;
    474637    }
    475  
    476 
    477   fRSigma    = TMath::Sqrt(rsigmasquare);
    478   fRSigmaVar = 0.25 * (sigmaSquareVar + pedRmsSquareVar) / rsigmasquare;
     638
     639  fRSigmaSquareVar = 4. * (sigmaSquareVar + pedRmsSquareVar);
    479640
    480641  return kTRUE;
    481642}
    482643
    483 //
    484 // Calculate the number of photo-electrons after the F-Factor method
    485 // Calculate the errors of the F-Factor method
     644// ------------------------------------------------------------------
     645//
     646// If fRSigmaSquare is smaller than 0 (i.e. has not yet been set),
     647// set kFFactorMethodValid to kFALSE and return kFALSE
     648//
     649// Calculate the number of photo-electrons with the F-Factor method:
     650// - Test bit IsHiGainSaturation() for the Mean Sum of FADC slices:
     651//   If yes, take fLoGainMean and fLoGainMeanVar
     652//   If no , take fHiGainMean and fHiGainMeanVar
     653//
     654// - Test bit IsHiGainSaturation() for the pedRMS:
     655//   If yes, take fLoGainPedRmsSquare and fLoGainPedRmsSquareVar
     656//   If no , take fPedRms and fPedRmsVar
     657//
     658// - Calculate the number of photo-electrons with the formula:
     659//   fPheFFactorMethod   = gkFFactor*gkFFactor * mean * mean  / fRSigmaSquare
     660//
     661// - Calculate the Variance on the photo-electrons with the formula:
     662//   fPheFFactorMethodVar =  (  4. * gkFFactorErr * gkFFactorErr / ( gkFFactor * gkFFactor )
     663//                            + 4. * fMeanVar                    / ( mean      * mean      )
     664//                            + 4. * fRSigmaVar                  / ( fRSigma   * fRSigma   )
     665//                            ) * fPheFFactor * fPheFFactor
     666// - if fPheFFactorMethod is less than fPheFFactorMethodLimit,
     667//   Set kFFactorMethodValid to kFALSE and return kFALSE
     668//   else: Set kFFactorMethodValid to kTRUE and return kTRUE
    486669//
    487670Bool_t MCalibrationChargePix::CalcFFactorMethod()
    488671{
    489672
    490   if (fRSigma < 0.)
     673  if (fRSigmaSquare < 0.)
    491674    {
    492675      SetFFactorMethodValid(kFALSE);
    … …  
    500683  // Square all variables in order to avoid applications of square root
    501684  //
    502   // First the relative error squares
    503   //
    504   const Float_t meanSquare        =     mean * mean;
    505   const Float_t meanSquareRelVar  = 4.* var/ meanSquare;
    506 
    507   const Float_t ffactorsquare       =    gkFFactor    * gkFFactor;
    508   const Float_t ffactorsquareRelVar = 4.*gkFFactorErr * gkFFactorErr / ffactorsquare;
    509 
    510   const Float_t rsigmaSquare        =     fRSigma    * fRSigma;
    511   const Float_t rsigmaSquareRelVar  = 4.* fRSigmaVar / rsigmaSquare;
     685  const Float_t meanSquare          =     mean * mean;
     686  const Float_t meanSquareRelVar    = 4.* var / meanSquare;
     687
     688  const Float_t ffactorsquare       =     gkFFactor    * gkFFactor;
     689  const Float_t ffactorsquareRelVar = 4.* gkFFactorErr * gkFFactorErr / ffactorsquare;
     690
     691  const Float_t rsigmaSquareRelVar  =     fRSigmaSquareVar / fRSigmaSquare;
    512692
    513693  //
    … …  
    515695  // (independent on Hi Gain or Lo Gain)
    516696  //
    517   fPheFFactorMethod = ffactorsquare * meanSquare / rsigmaSquare;
     697  fPheFFactorMethod = ffactorsquare * meanSquare / fRSigmaSquare;
    518698
    519699  if (fPheFFactorMethod < fPheFFactorMethodLimit)
    … …  
    534714
    535715
     716// ----------------------------------------------------------------------------
     717//
     718// - If fPed    is smaller than 0 (i.e. has not yet been set), return.
     719// - If fPedVar is smaller than 0 (i.e. has not yet been set), return.
     720//
     721// Calculate the electronic pedestal RMS with the formula:
     722//  - elec. pedestal = gkElectronicPedRms * sqrt(logainsamples)
     723//
     724// Calculate the LONS ped. RMS contribution in the high-gain
     725// from the high gain Pedestal RMS with the formula:
     726// - HiGain NSB square      = fPedRms * fPedRms - elec.ped.* elec.ped.
     727// - Var(HiGain NSB square) = fPedVar * fPedRms * fPedRms + 4.*elecPedRmsVar * elec.ped.* elec.ped.
     728//
     729// If PedRMS (LONS,lowgain) square is smaller than 0., set it to zero. (but not the error!)
     730//
     731// Convert the LONS ped. RMS contribution to the low-gain with the formula:
     732// - LoGain NSB square      = - HiGain NSB square / (fConversionHiLo*fConversionHiLo)
     733// - Var(LoGain NSB square) = ( Var(HiGain NSB square) / (HiGain NSB square * HiGain NSB square)
     734//                              + GetConversionHiLoRelVar()   )
     735//                            * LoGain NSB square * LoGain NSB square
     736//
     737// - Low Gain Ped RMS Square       = LoGain NSB square      + elec.ped. square
     738//   Var (Low Gain Ped RMS Square) = Var(LoGain NSB square) + Var(elec.ped. square)
     739//
    536740void MCalibrationChargePix::CalcLoGainPedestal(Float_t logainsamples)
    537741{
     742
     743  if (fPedRms < 0.)
     744    return;
     745
     746  if (fPedVar < 0.)
     747    return;
    538748
    539749  const Float_t elecPedRms     = gkElectronicPedRms    * TMath::Sqrt(logainsamples);
    … …  
    552762  const Float_t elecRmsSquareVar = 4.*elecPedRmsVar * elecRmsSquare;
    553763 
    554   Float_t nsbSquare             =  pedRmsSquare    - elecRmsSquare;
    555   Float_t nsbSquareRelVar       = (pedRmsSquareVar + elecRmsSquareVar)
    556                                  / (nsbSquare * nsbSquare) ;
    557  
    558   if (nsbSquare < 0.)
    559     nsbSquare = 0.;
     764  Float_t higainNsbSquare        =  pedRmsSquare    - elecRmsSquare;
     765  Float_t higainNsbSquareRelVar  = (pedRmsSquareVar + elecRmsSquareVar)
     766                                 / (higainNsbSquare * higainNsbSquare) ;
     767 
     768  if (higainNsbSquare < 0.)
     769    higainNsbSquare = 0.;
    560770 
    561771  //
    … …  
    563773  // add it quadratically to the electronic noise
    564774  //
    565   const Float_t conversionSquare        =    fConversionHiLo    * fConversionHiLo;
    566   const Float_t convertedNsbSquare      =    nsbSquare       / conversionSquare;
    567   const Float_t convertedNsbSquareVar   =    nsbSquareRelVar
    568                                             * convertedNsbSquare * convertedNsbSquare;
     775  const Float_t conversionSquare        =     fConversionHiLo    * fConversionHiLo;
     776  const Float_t conversionSquareRelVar  = 4.* GetConversionHiLoRelVar();
     777
     778  const Float_t logainNsbSquare      =      higainNsbSquare  / conversionSquare;
     779  const Float_t logainNsbSquareVar   =    ( higainNsbSquareRelVar + conversionSquareRelVar )
     780                                            * logainNsbSquare * logainNsbSquare;
    569781   
    570   pedRmsSquare     = convertedNsbSquare    + elecRmsSquare;
    571   pedRmsSquareVar  = convertedNsbSquareVar + elecRmsSquareVar;
    572  
    573   fLoGainPedRms    = TMath::Sqrt(pedRmsSquare);
    574   fLoGainPedRmsVar = 0.25 * pedRmsSquareVar /  pedRmsSquare;
    575  
     782  fLoGainPedRmsSquare    = logainNsbSquare    + elecRmsSquare;
     783  fLoGainPedRmsSquareVar = logainNsbSquareVar + elecRmsSquareVar;
    576784}
    577785 
     786const Float_t MCalibrationChargePix::GetConversionHiLoRelVar() const
     787{
     788  if (fConversionHiLo == 0.)
     789    return 0.;
     790
     791  return fConversionHiLoVar / (fConversionHiLo * fConversionHiLo);
     792}
    578793 
    579  

  • trunk/MagicSoft/Mars/mcalib/MHCalibrationCam.cc

    r3642 r3650  
    703703 
    704704  hist.Renorm();
     705
    705706  //
    706707  // 4) Check for oscillations
    … …  
    708709  hist.CreateFourierSpectrum();
    709710 
     711
    710712  if (!hist.IsFourierSpectrumOK())
    711713    bad.SetUncalibrated( osctyp );

  • trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeCam.cc

    r3639 r3650  
    216216{
    217217
    218   fCam = (MCalibrationCam*)pList->FindCreateObj("MCalibrationChargeCam");
     218  const Int_t npixels  = fGeom->GetNumPixels();
     219  const Int_t nsectors = fGeom->GetNumSectors();
     220  const Int_t nareas   = fGeom->GetNumAreas();
     221
     222  fCam = (MCalibrationCam*)pList->FindObject("MCalibrationChargeCam");
    219223  if (!fCam)
    220       return kFALSE;
     224    {
     225      fCam = (MCalibrationCam*)pList->FindCreateObj(AddSerialNumber("MCalibrationChargeCam"));
     226      if (!fCam)
     227        {
     228          gLog << err << "Cannot find nor create MCalibrationChargeCam ... abort." << endl;
     229          return kFALSE;
     230        }
     231      else
     232        {
     233          fCam->InitSize(npixels);
     234          fCam->InitAverageAreas(nareas);
     235          fCam->InitAverageSectors(nsectors);
     236        }
     237    }
    221238
    222239  MExtractedSignalCam *signal = (MExtractedSignalCam*)pList->FindObject("MExtractedSignalCam");
    … …  
    226243      return kFALSE;
    227244  }
    228 
    229   const Int_t npixels  = fGeom->GetNumPixels();
    230   const Int_t nsectors = fGeom->GetNumSectors();
    231   const Int_t nareas   = fGeom->GetNumAreas();
    232245
    233246  if (fHiGainArray->GetEntries()==0)

  • trunk/MagicSoft/Mars/mcalib/MHCalibrationRelTimeCam.cc

    r3642 r3650  
    153153{
    154154
    155   fCam = (MCalibrationCam*)pList->FindCreateObj("MCalibrationRelTimeCam");
     155  fCam = (MCalibrationCam*)pList->FindObject("MCalibrationRelTimeCam");
    156156  if (!fCam)
    157       return kFALSE;
    158 
     157    {
     158      fCam = (MCalibrationCam*)pList->FindCreateObj(AddSerialNumber("MCalibrationRelTimeCam"));
     159      if (!fCam)
     160        {
     161          gLog << err << "Cannot find nor create MCalibrationRelTimeCam ... abort." << endl;
     162          return kFALSE;
     163        }
     164      else
     165        {
     166          fCam->InitSize(npixels);
     167          fCam->InitAverageAreas(nareas);
     168          fCam->InitAverageSectors(nsectors);
     169        }
     170    }
     171 
    159172  MArrivalTimeCam *signal = (MArrivalTimeCam*)pList->FindObject("MArrivalTimeCam");
    160173  if (!signal)
    … …  
    376389{
    377390
    378   FitHiGainArrays((MCalibrationCam&)(*fCam),*fBadPixels,
     391  FitHiGainArrays((*fCam),*fBadPixels,
    379392                  MBadPixelsPix::kRelTimeNotFitted,
    380393                  MBadPixelsPix::kRelTimeOscillating);
    381394
    382   FitLoGainArrays((MCalibrationCam&)(*fCam),*fBadPixels,
     395  FitLoGainArrays((*fCam),*fBadPixels,
    383396                  MBadPixelsPix::kRelTimeNotFitted,
    384397                  MBadPixelsPix::kRelTimeOscillating);
Note: See TracChangeset for help on using the changeset viewer.