Changeset 2886


Ignore:
Timestamp:
01/22/04 20:54:16 (21 years ago)
Author:
moralejo
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • trunk/MagicSoft/Mars/Changelog

    r2885 r2886  
    44
    55                                                 -*-*- END OF LINE -*-*-
     6
     7 2004/01/22: Abelardo Moralejo
     8
     9   * manalysis/MMcCalibrationUpdate.[h,cc]
     10     - Now ratio of high to low gain is taken from MCalibrationCam if
     11       it existed previously in the parameter list, instead of being
     12       read again from the MMcFadcHeader. Removed Setter function for
     13       fADC2PhInner, no longer necessary. Fixed error regarding the
     14       pedestal conversion to photons (did not read conversion factor
     15       from preexisting MCalibrationCam object).
     16
     17   * mcalib/MMcCalibrationCalc.cc
     18     - Changed parameters of the histogram, and also the quantity being
     19       histogrammed. Check that input data come from a noiseless camera
     20       file before proceeding to do the calibration. Introduced lower
     21       size in cut for calibration. Now rhe calibration constant is not
     22       calculated from the mean of photons/ADC counts, but from the peak
     23       of the histogram.
     24
     25   * macros/starmc.C
     26     - Introduced new scheme. Now there are two loops over two different
     27       sets of files. First loop calculates the calibration constants,
     28       second one does the analysis. Introduced comments. Now the
     29       histogram used in the light calibration is written to the output
     30       file.
     31
    632 2004/01/22: Thomas Bretz
    733

  • trunk/MagicSoft/Mars/macros/starmc.C

    r2874 r2886  
    11/* ======================================================================== *\
    2 !
    3 ! *
    4 ! * This file is part of MARS, the MAGIC Analysis and Reconstruction
    5 ! * Software. It is distributed to you in the hope that it can be a useful
    6 ! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
    7 ! * It is distributed WITHOUT ANY WARRANTY.
    8 ! *
    9 ! * Permission to use, copy, modify and distribute this software and its
    10 ! * documentation for any purpose is hereby granted without fee,
    11 ! * provided that the above copyright notice appear in all copies and
    12 ! * that both that copyright notice and this permission notice appear
    13 ! * in supporting documentation. It is provided "as is" without express
    14 ! * or implied warranty.
    15 ! *
    16 !
    17 !
    18 !   Author(s): Abelardo Moralejo 1/2004 <mailto:moralejo@pd.infn.it>
    19 !              Thomas Bretz  5/2002 <mailto:tbretz@astro.uni-wuerzburg.de>
    20 !
    21 !   Copyright: MAGIC Software Development, 2000-2004
    22 !
    23 !
    24 \* ======================================================================== */
     2   !
     3   ! *
     4   ! * This file is part of MARS, the MAGIC Analysis and Reconstruction
     5   ! * Software. It is distributed to you in the hope that it can be a useful
     6   ! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
     7   ! * It is distributed WITHOUT ANY WARRANTY.
     8   ! *
     9   ! * Permission to use, copy, modify and distribute this software and its
     10   ! * documentation for any purpose is hereby granted without fee,
     11   ! * provided that the above copyright notice appear in all copies and
     12   ! * that both that copyright notice and this permission notice appear
     13   ! * in supporting documentation. It is provided "as is" without express
     14   ! * or implied warranty.
     15   ! *
     16   !
     17   !
     18   !   Author(s): Abelardo Moralejo 1/2004 <mailto:moralejo@pd.infn.it>
     19   !              Thomas Bretz  5/2002 <mailto:tbretz@astro.uni-wuerzburg.de>
     20   !
     21   !   Copyright: MAGIC Software Development, 2000-2004
     22   !
     23   !
     24   \* ======================================================================== */
    2525
    2626/////////////////////////////////////////////////////////////////////////////
    … …  
    3737void starmc()
    3838{
    39     //
    40     // This is a demonstration program which calculates the image
    41     // parameters from Magic Monte Carlo files (output of camera).
    42 
    43     //
    44     // Create a empty Parameter List and an empty Task List
    45     // The tasklist is identified in the eventloop by its name
    46     //
    47     MParList  plist;
    48 
    49     MTaskList tlist;
    50 
    51     plist.AddToList(&tlist);
    52 
    53     MSrcPosCam src;
    54     src.SetReadyToSave();
    55 
    56     plist.AddToList(&src);
    57 
    58     //
    59     // Now setup the tasks and tasklist:
    60     // ---------------------------------
    61     //
    62     MReadMarsFile read("Events");
    63 
    64     // ------------- user change -----------------
    65 
    66     read.AddFile("Gamma_zbin*.root");
    67 
    68     read.DisableAutoScheme();
    69 
    70     MGeomApply geom; // Reads in geometry from MC file and sets the right sizes for
    71                      // several parameter containers.
    72 
    73     MMcPedestalCopy   pcopy;
    74     // Copies pedestal data from the MC file run fadc header to the MPedestalCam container.
    75 
    76     MExtractSignal    sigextract;
    77     // Define ADC slices to be integrated in high and low gain:
    78     sigextract.SetRange(0, 5, 0, 5);
    79 
    80     MMcCalibrationUpdate  mccalibupdate;
    81 
    82     //
    83     // Now introduce conversion factor from ADC counts to photons before camera for
    84     // inner pixels. The corresponding value for outer pixels is then calculated
    85     // automatically. Bear in mind that the conversion factor depend both on the
    86     // parameters used in the camera simulation as well as on the number of
    87     // integrated FADC slices. In the future it should be calculated in a previous
    88     // event loop, either from the same MC file or from a MC calibration file
    89     // written on purpose.
    90     // (FIXME: the conversion must be calculated automatically from the analyzed file)
    91     //
    92     mccalibupdate.SetADC2PhInner(1.2586);
     39  //
     40  // This is a demonstration program which calculates the image
     41  // parameters from Magic Monte Carlo files (output of camera).
    9342
    9443
    95     MCalibrate calib; // Transforms signals from ADC counts into photons.
     44  // ------------- user change -----------------
    9645
    97     //    MBlindPixelCalc   blind;
    98     //    blind.SetUseInterpolation();
     46  TString* CalibrationFilename;
    9947
    100     MImgCleanStd      clean(4.1,3.4); // Applies tail cuts to image.
     48  //
     49  // Comment out next line to disable calibration. In that case the units of the
     50  // MHillas.fSize parameter will be ADC counts (rather, equivalent ADC counts in
     51  // inner pixels, since we correct for the possible differences in gain of outer
     52  // pixels)
     53  //
     54  CalibrationFilename = new TString("../../../nonoise/gammas/Gamma_zbin0_0*.root");
     55  // File to be used in the calibration (must be a camera file without added noise)
     56
     57  Char_t* AnalysisFilename = "Proton_zbin*.root";  // File to be analyzed
     58  Char_t* OutFilename      = "star_mc.root";       // Output file name
     59
     60  Float_t CleanLev[2] = {4., 3.}; // Tail cuts for image analysis
     61
     62  Int_t BinsHigh[2] = {0, 5}; // First and last FADC bin of the range to be integrated,
     63  Int_t BinsLow[2]  = {0, 5}; // for high and low gain respectively.
     64
     65  // -------------------------------------------
     66
     67  //
     68  // Create a empty Parameter List and an empty Task List
     69  // The tasklist is identified in the eventloop by its name
     70  //
     71  MParList  plist;
     72
     73  MTaskList tlist;
     74
     75  plist.AddToList(&tlist);
     76
     77  MSrcPosCam src;
     78  src.SetReadyToSave();
     79
     80  plist.AddToList(&src);
     81
     82  //
     83  // Now setup the tasks and tasklist:
     84  // ---------------------------------
     85  //
     86  MReadMarsFile read("Events");
     87
     88  if (CalibrationFilename)
     89    read.AddFile(CalibrationFilename->Data());
     90
     91  read.DisableAutoScheme();
     92
     93  MGeomApply geom; // Reads in geometry from MC file and sets the right sizes for
     94  // several parameter containers.
     95
     96  MMcPedestalCopy   pcopy;
     97  // Copies pedestal data from the MC file run fadc header to the MPedestalCam container.
     98
     99  MExtractSignal    sigextract;
     100  // Define ADC slices to be integrated in high and low gain:
     101  sigextract.SetRange(BinsHigh[0], BinsHigh[1], BinsLow[0], BinsLow[1]);
     102
     103  MMcCalibrationUpdate  mccalibupdate;
     104
     105  MCalibrate calib; // Transforms signals from ADC counts into photons.
     106
     107  //    MBlindPixelCalc   blind;
     108  //    blind.SetUseInterpolation();
     109
     110  MImgCleanStd      clean(CleanLev[0], CleanLev[1]); // Applies tail cuts to image.
     111
     112  MHillasCalc       hcalc; // Calculates Hillas parameters not dependent on source position.
     113  MHillasSrcCalc    scalc; // Calculates source-dependent Hillas parameters
     114
     115  MMcCalibrationCalc mccalibcalc;
     116
     117  tlist.AddToList(&read);
     118  tlist.AddToList(&geom);
     119  tlist.AddToList(&pcopy);
     120
     121  tlist.AddToList(&sigextract);
     122  tlist.AddToList(&mccalibupdate);
     123  tlist.AddToList(&calib);
     124  tlist.AddToList(&clean);
     125  //    tlist.AddToList(&blind);
     126  tlist.AddToList(&hcalc);
     127
     128  tlist.AddToList(&mccalibcalc);
     129
     130  //
     131  // Open output file:
     132  //
     133  MWriteRootFile write(OutFilename); // Writes output
     134  write.AddContainer("MRawRunHeader", "RunHeaders");
     135  write.AddContainer("MMcRunHeader",  "RunHeaders");
     136  write.AddContainer("MSrcPosCam",    "RunHeaders");
     137  write.AddContainer("MMcEvt",        "Events");
     138  write.AddContainer("MHillas",       "Events");
     139  write.AddContainer("MHillasExt",    "Events");
     140  write.AddContainer("MHillasSrc",    "Events");
     141  write.AddContainer("MNewImagePar",  "Events");
     142
     143  //
     144  // First loop: Calibration loop
     145  //
     146
     147  MProgressBar bar;
     148  bar.SetWindowName("Calibrating...");
     149
     150  MEvtLoop evtloop;
     151  evtloop.SetProgressBar(&bar);
     152  evtloop.SetParList(&plist);
     153
     154  if (CalibrationFilename)
     155    {
     156      if (!evtloop.Eventloop())
     157        return;
     158      mccalibcalc->GetHist()->Write();
     159    }
     160
     161  //
     162  // Second loop: analysis loop
     163  //
     164
     165  //
     166  // Change the read task by another one which reads the file we want to analyze:
     167  //
     168
     169  MReadMarsFile read2("Events");
     170  read2.AddFile(AnalysisFilename);
     171  read2.DisableAutoScheme();
     172  tlist.AddToListBefore(&read2, &read, "All");
     173  tlist.RemoveFromList(&read);
     174
     175  bar.SetWindowName("Analyzing...");
     176
     177  tlist.RemoveFromList(&mccalibcalc); // Removes calibration task from list.
     178
     179  tlist.AddToList(&scalc);            // Calculates Source-dependent Hillas parameters
     180
     181  tlist.AddToList(&write);            // Add task to write output.
     182
     183  if (!evtloop.Eventloop())
     184    return;
    101185
    102186
    103     MHillasCalc       hcalc; // Calculates Hillas parameters not dependent on source position.
    104     MHillasSrcCalc    scalc; // Calculates source-dependent Hillas parameters
    105                              // (!!Preliminary!! Will be removed later!)
    106 
    107 
    108     // ------------- user change -----------------
    109 
    110     MWriteRootFile write("star_mc.root"); // Writes output
    111 
    112     write.AddContainer("MRawRunHeader", "RunHeaders");
    113     write.AddContainer("MMcRunHeader",  "RunHeaders");
    114     write.AddContainer("MSrcPosCam",    "RunHeaders");
    115     write.AddContainer("MMcEvt",        "Events");
    116     write.AddContainer("MHillas",       "Events");
    117     write.AddContainer("MHillasExt",    "Events");
    118     write.AddContainer("MHillasSrc",    "Events");
    119     write.AddContainer("MNewImagePar",  "Events");
    120 
    121 
    122     tlist.AddToList(&read);
    123     tlist.AddToList(&geom);
    124     tlist.AddToList(&pcopy);
    125 
    126     tlist.AddToList(&sigextract);
    127     tlist.AddToList(&mccalibupdate);
    128     tlist.AddToList(&calib);
    129     tlist.AddToList(&clean);
    130     //    tlist.AddToList(&blind);
    131     tlist.AddToList(&hcalc);
    132     tlist.AddToList(&scalc);
    133     tlist.AddToList(&write);
    134 
    135     //
    136     // Create and set up the eventloop
    137     //
    138     MProgressBar bar;
    139 
    140     MEvtLoop evtloop;
    141     evtloop.SetProgressBar(&bar);
    142     evtloop.SetParList(&plist);
    143 
    144     //
    145     // Execute your analysis
    146      //
    147     if (!evtloop.Eventloop())
    148         return;
    149 
    150     tlist.PrintStatistics();
     187  tlist.PrintStatistics();
    151188}

  • trunk/MagicSoft/Mars/manalysis/MMcCalibrationUpdate.cc

    r2876 r2886  
    3737//   MMcFadcHeader
    3838//   MRawRunHeader
     39//  [MCalibrationCam] (if it existed previously)
    3940//
    4041//  Output Containers:
    41 //   MCalibrationCam
    4242//   MPedPhotCam
     43//  [MCalibrationCam] (if it did not exist previously)
    4344//
    4445/////////////////////////////////////////////////////////////////////////////
    … …  
    9596}
    9697
    97 
    98 //---------------------------------------------------------------------------
    99 //
    100 // Set ADC to photon conversion factors.
    101 //
    102 void MMcCalibrationUpdate::SetADC2PhInner(Float_t x)
    103 {
    104   fADC2PhInner = x;
    105 
    106   return;
    107 }
    108 
    10998// --------------------------------------------------------------------------
    11099//
    … …  
    123112  if ( !fCalCam )
    124113    {
    125       *fLog << warn << dbginf << AddSerialNumber("MCalibrationCam") << " does not exist... Creating." << endl;
     114      *fLog << inf << dbginf << AddSerialNumber("MCalibrationCam") << " does not exist... Creating." << endl;
    126115
    127116      fCalCam = (MCalibrationCam*) pList->FindCreateObj(AddSerialNumber("MCalibrationCam"));
    … …  
    135124    {
    136125      fFillCalibrationCam = kFALSE;
    137       *fLog << warn << dbginf << AddSerialNumber("MCalibrationCam") << " already exists... " << endl;
     126      *fLog << inf << AddSerialNumber("MCalibrationCam") << " already exists... " << endl;
    138127    }
    139128
    … …  
    189178
    190179    //
    191     // Set the ADC to photons conversion factor for outer pixels:
    192     //
    193     fADC2PhOuter = fADC2PhInner *
    194       (fHeaderFadc->GetAmplitud() / fHeaderFadc->GetAmplitudOuter());
    195 
    196     //
    197     // Set the conversion factor between high and low gain:
    198     //
    199     fConversionHiLo = fHeaderFadc->GetLow2HighGain();
    200 
    201     //
    202180    // If MCalibrationCam already existed in the parameter list before
    203181    // MMcCalibrationUpdate::PreProcess was executed (from a
    … …  
    208186      return kTRUE;
    209187
     188    //
     189    // Set the ADC to photons conversion factor for outer pixels:
     190    //
     191    fADC2PhOuter = fADC2PhInner *
     192      (fHeaderFadc->GetAmplitud() / fHeaderFadc->GetAmplitudOuter());
     193
     194    //
     195    // Set the conversion factor between high and low gain:
     196    //
     197    fConversionHiLo = fHeaderFadc->GetLow2HighGain();
    210198
    211199    const int num = fCalCam->GetSize();
    … …  
    264252        // counts for the RMS per slice:
    265253        //
    266 
    267254
    268255        const Float_t pedestrms  = sigpix.IsLoGainUsed()?
    … …  
    279266        MPedPhotPix &pedpix = (*fPedPhotCam)[i];
    280267
    281         Float_t adc2phot = (fGeom->GetPixRatio(i) < fGeom->GetPixRatio(0))?
    282           fADC2PhOuter : fADC2PhInner;
     268        MCalibrationPix &calpix = (*fCalCam)[i];
     269        Float_t adc2phot = calpix.GetMeanConversionBlindPixelMethod();
     270        Float_t hi2lo    = calpix.GetConversionHiLo();
    283271
    284272        if (sigpix.IsLoGainUsed())
    285           pedpix.Set(adc2phot*fConversionHiLo*pedestmean,
    286                      adc2phot*fConversionHiLo*pedestrms);
     273          pedpix.Set(adc2phot*hi2lo*pedestmean,
     274                     adc2phot*hi2lo*pedestrms);
    287275        else
    288276          pedpix.Set(adc2phot*pedestmean, adc2phot*pedestrms);

  • trunk/MagicSoft/Mars/manalysis/MMcCalibrationUpdate.h

    r2876 r2886  
    3535    MMcCalibrationUpdate(const char *name=NULL, const char *title=NULL);
    3636
    37     void SetADC2PhInner(Float_t x);
    38 
    3937    ClassDef(MMcCalibrationUpdate, 0)   // Task which obtains, for MC files, the pedestal mean and rms, and the calibration factor from ADC counts to photons.
    4038};

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

    r2885 r2886  
    7171  fEvents = 0;
    7272
    73   fHistRatio = new TH1F("HistRatio", "Ratio fPassPhotCone/fSize", 1000., 0., 10.);
     73  fHistRatio = new TH1F("HistRatio", "log10(fPassPhotCone/fSize)", 1500., -3., 3.);
     74  fHistRatio->GetXaxis()->SetTitle("log_{10}(fPassPhotCone / fSize) (in photons/ADC count)");
    7475}
    7576
    … …  
    165166    }
    166167
     168  for (UInt_t ipix = 0; ipix < fGeom->GetNumPixels(); ipix++)
     169    {
     170      if (fHeaderFadc->GetPedestalRmsHigh(ipix) > 0 ||
     171          fHeaderFadc->GetPedestalRmsLow(ipix) > 0 )
     172        {
     173          *fLog << err << endl << endl << dbginf << "You are trying to calibrate the data using a Camera file produced with added noise. Please use a noiseless file for calibration. Aborting..." << endl << endl;
     174          return kFALSE;
     175        }
     176    }
     177
    167178  //
    168179  // FIXME: Check that the relevant parameters in the FADC header are the
    169   // same for all the analyzed data!
     180  // same for all the analyzed data, both Calibration and Analyzed data!
    170181  //
    171182
    … …  
    181192{
    182193
     194  //
     195  // Exclude events with some saturated pixel
     196  //
    183197  if ( fNew->GetNumSaturatedPixels() > 0 )
    184198    return kTRUE;
    185199
     200  //
     201  // Exclude events with low Size (larger fluctuations)
     202  // FIXME? The present cut (1000 "inner-pixel-counts") is somehow arbitrary.
     203  // Might it be optimized?
     204  //
     205  if ( fHillas->GetSize() < 1000 )
     206    return kTRUE;
     207
    186208  fADC2Phot += fMcEvt->GetPassPhotCone() / fHillas->GetSize();
    187209  fEvents ++;
    188210
    189   fHistRatio->Fill(fMcEvt->GetPassPhotCone()/fHillas->GetSize());
     211  fHistRatio->Fill(log10(fMcEvt->GetPassPhotCone()/fHillas->GetSize()));
    190212
    191213  return kTRUE;
    … …  
    207229    }
    208230
     231  //
     232  // For the calibration we no longer use the mean, but thepeak of the distribution:
     233  //
     234
     235  Float_t summax = 0.;
     236  Int_t mode = 0;
     237  Int_t reach = 2;
     238  for (Int_t ibin = 1+reach; ibin <= fHistRatio->GetNbinsX()-reach; ibin++)
     239    {
     240      Float_t sum = 0;
     241      for(Int_t k = ibin-reach; k <= ibin+reach; k++)
     242        sum += fHistRatio->GetBinContent(k);
     243      if (sum > summax)
     244        {
     245          summax = sum;
     246          mode = ibin;
     247        }
     248    }
     249
     250  fADC2Phot = pow(10., fHistRatio->GetXaxis()->GetBinCenter(mode));
     251
    209252  const int num = fCalCam->GetSize();
    210253 

  • trunk/MagicSoft/Mars/mcalib/MMcCalibrationCalc.h

    r2885 r2886  
    2828    Long_t  fEvents;
    2929
    30     TH1F*   fHistRatio;
     30    TH1F*   fHistRatio; // Histogram for monitoring the calibration.
    3131
    3232    Bool_t CheckRunType(MParList *pList) const;
    … …  
    3939    MMcCalibrationCalc(const char *name=NULL, const char *title=NULL);
    4040
     41    TH1F*   GetHist() { return fHistRatio; }
     42
    4143    ClassDef(MMcCalibrationCalc, 0)   // Task which obtains, for MC files, the calibration factor from ADC counts to photons.
    4244};
Note: See TracChangeset for help on using the changeset viewer.