Ignore:
Timestamp:
02/12/05 15:56:23 (20 years ago)
Author:
moralejo
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars/macros
Files:
3 edited

Legend:

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

    r6357 r6403  
    2929//
    3030//  This macro converts raw MC data into calibrated data (photons per pixel)
     31//  It optionally divides the output into a train and a test sample
    3132//
    3233/////////////////////////////////////////////////////////////////////////////
     
    3738{
    3839  //
    39   // This is a demonstration program which reads in MC camera files
    40   // and produces and output with calibrated events (signal in photons
    41   // for all pixels, in MCerPhotEvt containers).
     40  // This macro reads in MC camera files and produces and output with
     41  // calibrated events (signal in photons for all pixels, in MCerPhotEvt
     42  // containers).
    4243  //
    4344
     
    4748  CalibrationFilename = new TString("/users/emc/moralejo/mcdata/Period021_0.73_mirror/gammas_nonoise/Gamma_*root");  // File to be used for the calibration (must be a camera file without added noise)
    4849
    49   Char_t* AnalysisFilename = "Gamma_*.root";  // File to be analyzed
    50 
    51   Char_t* OutFilename      = "calibrated_gamma.root";  // Output file name
    52 
    53 
    54   // (other extraction methods can be used)
     50  Char_t* AnalysisFilename = "Gamma_zbin*w0.root";  // File to be analyzed
     51
     52
     53  TString* OutFilename1;
     54  TString* OutFilename2;
     55
     56  // Output file names
     57  OutFilename1 = new TString("calibrated_gamma_train.root");
     58  OutFilename2 = new TString("calibrated_gamma_test.root");
     59
     60  // To get only one output file, just comment out the second
     61  // one of the above lines
     62
     63
     64  //
     65  //  Set signal extractor
     66  //
    5567  //    MExtractFixedWindowPeakSearch sigextract;
    5668  //    sigextract.SetWindows(6, 6, 4);
    57 
     69  //
    5870  MExtractTimeAndChargeDigitalFilter sigextract;
    5971  sigextract.SetNameWeightsFile("/users/emc/moralejo/Mars/msignal/MC_weights.dat");
     
    100112
    101113  MPointingPosCalc pointcalc;
    102   // Creates MPointingPos object and fill it with the telescope orientation
    103   // information taken from MMcEvt.
     114  // Creates MPointingPos object and fills it with the telescope
     115  // orientation information taken from MMcEvt.
    104116
    105117  MCalibrateData calib;
    106   // MCalibrateData transforms signals from ADC counts into photons. In the first
    107   // loop it applies a "dummy" calibration supplied by MMcCalibrationUpdate, just
    108   // to equalize inner and outer pixels. At the end of the first loop, in the
    109   // PostProcess of MMcCalibrationCalc (see below) the true calibration constants
     118  //
     119  // MCalibrateData transforms signals from ADC counts into photons or phe-
     120  // (this can be selected anove). In the first loop it applies a "dummy"
     121  // calibration supplied by MMcCalibrationUpdate, just to equalize inner
     122  // and outer pixels, and calculates SIZE in "equivalent number of inner
     123  // ADC counts". At the end of the first loop, in the PostProcess of
     124  // MMcCalibrationCalc (see below) the true calibration constants
    110125  // are calculated.
    111 
     126  //
    112127  calib.SetCalibrationMode(MCalibrateData::kFfactor);
    113128  // Do not change the CalibrationMode above for MC...!
     
    124139  //
    125140
    126   MHillasCalc hcalc; // Calculates Hillas parameters not dependent on source position.
     141  MHillasCalc hcalc;
    127142  hcalc.Disable(MHillasCalc::kCalcHillasSrc);
     143  // Calculates Hillas parameters not dependent on source position.
    128144
    129145  MMcCalibrationCalc mccalibcalc;
    130146  // Calculates calibration constants to convert from ADC counts to photons.
    131  
    132147
    133148  tlist.AddToList(&read);
    134149  tlist.AddToList(&geom);
     150
     151  MF notrigger("MMcTrig.fNumFirstLevel<1");
     152  MContinue skipnotrig(&notrigger);
     153  tlist.AddToList(&skipnotrig);
     154
    135155  tlist.AddToList(&pcopy);
    136156  tlist.AddToList(&pointcalc);
     
    144164
    145165  //
    146   // Open output file:
    147   //
    148   MWriteRootFile write(OutFilename); // Writes output
    149   write.AddContainer("MGeomCam",            "RunHeaders");
    150   write.AddContainer("MMcConfigRunHeader",  "RunHeaders");
    151   write.AddContainer("MMcCorsikaRunHeader", "RunHeaders");
    152   write.AddContainer("MMcFadcHeader",       "RunHeaders");
    153   write.AddContainer("MMcRunHeader",        "RunHeaders");
    154   write.AddContainer("MMcTrigHeader",       "RunHeaders");
    155   write.AddContainer("MRawRunHeader",       "RunHeaders");
    156 
    157 
    158   write.AddContainer("MMcEvt",        "Events");
    159   write.AddContainer("MMcTrig",       "Events");
    160   write.AddContainer("MPointingPos",  "Events");
    161   write.AddContainer("MRawEvtHeader", "Events");
    162   write.AddContainer("MCerPhotEvt",   "Events");
    163   write.AddContainer("MPedPhotCam",   "Events");
     166  // Open output files:
     167  //
     168  MWriteRootFile write1(OutFilename1->Data()); // Writes output
     169
     170  MWriteRootFile write1_b(OutFilename1->Data());
     171  write1_b.AddContainer("MMcEvtBasic", "OriginalMC", kFALSE);
     172  // Add the MMcEvtBasic container only in case it exists in
     173  // the input camera files
     174
     175  write1.AddContainer("MGeomCam",            "RunHeaders");
     176  write1.AddContainer("MMcConfigRunHeader",  "RunHeaders");
     177  write1.AddContainer("MMcCorsikaRunHeader", "RunHeaders");
     178  write1.AddContainer("MMcFadcHeader",       "RunHeaders");
     179  write1.AddContainer("MMcRunHeader",        "RunHeaders");
     180  write1.AddContainer("MMcTrigHeader",       "RunHeaders");
     181  write1.AddContainer("MRawRunHeader",       "RunHeaders");
     182
     183  write1.AddContainer("MMcEvt",        "Events");
     184  write1.AddContainer("MMcTrig",       "Events");
     185  write1.AddContainer("MPointingPos",  "Events");
     186  write1.AddContainer("MRawEvtHeader", "Events");
     187  write1.AddContainer("MCerPhotEvt",   "Events");
     188  write1.AddContainer("MPedPhotCam",   "Events");
     189
     190  if (OutFilename2)
     191    {
     192      MWriteRootFile write2(OutFilename2->Data()); // Writes output   
     193
     194      MWriteRootFile write2_b(OutFilename2->Data());
     195      write2_b.AddContainer("MMcEvtBasic", "OriginalMC", kFALSE);
     196      // Add the MMcEvtBasic container only in case it exists in
     197      // the input camera files
     198
     199      write2.AddContainer("MGeomCam",            "RunHeaders");
     200      write2.AddContainer("MMcConfigRunHeader",  "RunHeaders");
     201      write2.AddContainer("MMcCorsikaRunHeader", "RunHeaders");
     202      write2.AddContainer("MMcFadcHeader",       "RunHeaders");
     203      write2.AddContainer("MMcRunHeader",        "RunHeaders");
     204      write2.AddContainer("MMcTrigHeader",       "RunHeaders");
     205      write2.AddContainer("MRawRunHeader",       "RunHeaders");
     206
     207      write2.AddContainer("MMcEvt",        "Events");
     208      write2.AddContainer("MMcTrig",       "Events");
     209      write2.AddContainer("MPointingPos",  "Events");
     210      write2.AddContainer("MRawEvtHeader", "Events");
     211      write2.AddContainer("MCerPhotEvt",   "Events");
     212      write2.AddContainer("MPedPhotCam",   "Events");
     213
     214      //
     215      // Divide output in train and test samples, using the event number
     216      // (odd/even) to achieve otherwise unbiased event samples:
     217      //
     218      MF filter1("{MMcEvt.fEvtNumber%2}>0.5");
     219      MF filter2("{MMcEvt.fEvtNumber%2}<0.5");
     220
     221      write1.SetFilter(&filter1);
     222      write2.SetFilter(&filter2);
     223
     224      write1_b.SetFilter(&filter1);
     225      write2_b.SetFilter(&filter2);
     226    }
    164227
    165228  //
     
    180243      mccalibcalc->GetHistADC2PhotEl()->Write();
    181244      mccalibcalc->GetHistPhot2PhotEl()->Write();
    182       // Writes out the histograms used for calibration.
     245      // Writes out the histograms used for calibration. In case of
     246      // aslit output in train and test file, this is written only
     247      // in the test file for now.
    183248    }
    184249
     
    198263  tlist.RemoveFromList(&read);
    199264
     265  // Now add tasks to write MMcEvtBasic to the OriginalMC tree:
     266  tlist.AddToListBefore(&write1_b, &skipnotrig, "All");
     267  if (OutFilename2)
     268    tlist.AddToListBefore(&write2_b, &skipnotrig, "All");
     269
     270
     271  if (OutFilename2) // Add filters to split output in train and test
     272    {
     273      tlist.AddToListBefore(&filter1, &write1_b, "All");
     274      tlist.AddToListBefore(&filter2, &write1_b, "All");
     275    }
     276
    200277  bar.SetWindowName("Writing...");
    201278
     
    204281  tlist.RemoveFromList(&mccalibcalc);
    205282
    206   tlist.AddToList(&write);            // Add task to write output.
     283  tlist.AddToList(&write1);
     284  tlist.AddToList(&write2);
     285  // Add tasks to write the Events and RunHeaders trees to output.
    207286
    208287  if (!evtloop.Eventloop())
  • trunk/MagicSoft/Mars/macros/starmc.C

    r6342 r6403  
    7474  Float_t accepted_fraction = 1.;
    7575
    76   Float_t CleanLev[2] = {5.75, 3.84};
     76  Float_t CleanLev[2] = {5., 4.};
    7777  // User change: tail cuts for image analysis
    7878
     
    175175  //
    176176
    177   MWriteRootFile write1(OutFilename1.Data()); // Writes output1
     177  MWriteRootFile write1(OutFilename1->Data()); // Writes output1
    178178  write1.AddContainer("MRawRunHeader", "RunHeaders");
    179179  write1.AddContainer("MMcRunHeader",  "RunHeaders");
     
    197197  if (OutFilename2)
    198198    {
    199       MWriteRootFile write2(OutFilename2.Data()); // Writes output2
     199      MWriteRootFile write2(OutFilename2->Data()); // Writes output2
    200200      write2.AddContainer("MRawRunHeader", "RunHeaders");
    201201      write2.AddContainer("MMcRunHeader",  "RunHeaders");
  • trunk/MagicSoft/Mars/macros/starmc2.C

    r6357 r6403  
    3838void starmc2()
    3939{
    40   Char_t* AnalysisFilename = "calibrated_gamma.root"; // File to be analyzed
     40  Char_t* AnalysisFilename = "calibrated_gamma_train.root"; // File to be analyzed
     41  //  Char_t* AnalysisFilename = "calibrated_gamma_test.root"; // File to be analyzed
    4142
    42   TString* OutFilename1;
    43   TString* OutFilename2;
     43  TString* OutFilename;
    4444
    45   // Change output file names as desired. If you want only one output, comment
    46   // out the initialization of OutFilename2:
    47   OutFilename1 = new TString("star_gamma_train.root");   // Output file name 1 (test)
    48   OutFilename2 = new TString("star_gamma_test.root");    // Output file name 2 (train)
     45  // Output file name
     46  OutFilename = new TString("star_gamma_train.root");   // Output file name
     47  //  OutFilename = new TString("star_gamma_test.root");   // Output file name
    4948
    5049  MImgCleanStd      clean(4.5, 3.); // Applies tail cuts to image.
     
    6766  //
    6867  // FOR WOBBLE MODE!! Set source position on camera here.
    69   //  src.SetX(120.);  // units: mm
     68  //  src.SetX(120.); 
     69  // units: mm. This 120 mm correspond to MC wobble mode w+ for zbin>0
     70  //
    7071
    7172  src.SetReadyToSave();
     
    8081
    8182  read.AddFile(AnalysisFilename);
    82 
    8383  read.DisableAutoScheme();
    8484
     
    9494  // Open output file(s):
    9595  //
    96   MWriteRootFile write1(OutFilename1->Data()); // Writes output   
     96  MWriteRootFile write(OutFilename->Data()); // Writes output   
    9797  //
    98   // Add MC containers only if they exist. In this way you can also run on real calibrated data.
     98  // Add MC containers only if they exist.
     99  // In this way you can also run on real calibrated data.
    99100  //
    100   write1.AddContainer("MRawRunHeader", "RunHeaders");
    101   write1.AddContainer("MMcRunHeader",  "RunHeaders", kFALSE);
    102   write1.AddContainer("MGeomCam",      "RunHeaders", kFALSE);
    103   write1.AddContainer("MMcConfigRunHeader",  "RunHeaders", kFALSE);
    104   write1.AddContainer("MMcCorsikaRunHeader", "RunHeaders", kFALSE);
    105   write1.AddContainer("MMcFadcHeader",  "RunHeaders", kFALSE);
    106   write1.AddContainer("MMcTrigHeader",  "RunHeaders", kFALSE);
     101  write.AddContainer("MRawRunHeader", "RunHeaders");
     102  write.AddContainer("MMcRunHeader",  "RunHeaders", kFALSE);
     103  write.AddContainer("MGeomCam",      "RunHeaders", kFALSE);
     104  write.AddContainer("MMcConfigRunHeader",  "RunHeaders", kFALSE);
     105  write.AddContainer("MMcCorsikaRunHeader", "RunHeaders", kFALSE);
     106  write.AddContainer("MMcFadcHeader",  "RunHeaders", kFALSE);
     107  write.AddContainer("MMcTrigHeader",  "RunHeaders", kFALSE);
     108       
     109  write.AddContainer("MMcEvt",        "Events", kFALSE);
     110  write.AddContainer("MPointingPos",  "Events", kFALSE);
     111  write.AddContainer("MMcTrig",       "Events", kFALSE);
     112  write.AddContainer("MSrcPosCam",    "Events", kFALSE);
     113  write.AddContainer("MRawEvtHeader", "Events");
     114  write.AddContainer("MHillas",       "Events");
     115  write.AddContainer("MHillasExt",    "Events");
     116  write.AddContainer("MHillasSrc",    "Events");
     117  write.AddContainer("MImagePar",     "Events");
     118  write.AddContainer("MNewImagePar",  "Events");
     119  write.AddContainer("MConcentration","Events");
    107120
    108   write1.AddContainer("MMcEvt",        "Events", kFALSE);
    109   write1.AddContainer("MPointingPos",  "Events", kFALSE);
    110   write1.AddContainer("MMcTrig",       "Events", kFALSE);
    111   write1.AddContainer("MSrcPosCam",    "Events", kFALSE);
    112   write1.AddContainer("MRawEvtHeader", "Events");
    113   write1.AddContainer("MHillas",       "Events");
    114   write1.AddContainer("MHillasExt",    "Events");
    115   write1.AddContainer("MHillasSrc",    "Events");
    116   write1.AddContainer("MImagePar",     "Events");
    117   write1.AddContainer("MNewImagePar",  "Events");
    118   write1.AddContainer("MConcentration","Events");
    119 
    120   if (OutFilename2) // Second output file, in case we want a split output
    121     {
    122       MWriteRootFile write2(OutFilename2->Data()); // Writes output
    123       write2.AddContainer("MRawRunHeader", "RunHeaders");
    124       write2.AddContainer("MMcRunHeader",  "RunHeaders", kFALSE);
    125       write2.AddContainer("MGeomCam",      "RunHeaders", kFALSE);
    126       write2.AddContainer("MMcConfigRunHeader",  "RunHeaders", kFALSE);
    127       write2.AddContainer("MMcCorsikaRunHeader", "RunHeaders", kFALSE);
    128       write2.AddContainer("MMcFadcHeader",  "RunHeaders", kFALSE);
    129       write2.AddContainer("MMcTrigHeader",  "RunHeaders", kFALSE);
    130 
    131       write2.AddContainer("MMcEvt",        "Events", kFALSE);
    132       write2.AddContainer("MPointingPos",  "Events", kFALSE);
    133       write2.AddContainer("MMcTrig",       "Events", kFALSE);
    134       write2.AddContainer("MSrcPosCam",    "Events", kFALSE);
    135       write2.AddContainer("MRawEvtHeader", "Events");
    136       write2.AddContainer("MHillas",       "Events");
    137       write2.AddContainer("MHillasExt",    "Events");
    138       write2.AddContainer("MHillasSrc",    "Events");
    139       write2.AddContainer("MImagePar",     "Events");
    140       write2.AddContainer("MNewImagePar",  "Events");
    141       write2.AddContainer("MConcentration","Events");
    142 
    143       //
    144       // Divide output in train and test samples, using the event number
    145       // (odd/even) to achieve otherwise unbiased event samples:
    146       //
    147      
    148       MF filter1("{MMcEvt.fEvtNumber%2}>0.5");
    149       MF filter2("{MMcEvt.fEvtNumber%2}<0.5");
    150 
    151       write1.SetFilter(&filter1);
    152       write2.SetFilter(&filter2);
    153 
    154       tlist.AddToList(&filter1);
    155       tlist.AddToList(&filter2);
    156     }
    157 
    158 
    159   tlist.AddToList(&write1);
    160 
    161   if (OutFilename2)
    162     tlist.AddToList(&write2);   
    163  
     121  tlist.AddToList(&write);
    164122
    165123  //
     
    178136  tlist.PrintStatistics();
    179137
     138  ////////////////////////////////////////////////////////////////////
     139  //
     140  // Second event loop: simply copy the tree MMcEvtBasic in the tree
     141  // "OriginalMC" from the input file to the output file:
     142
     143  MParList  plist2;
     144  MTaskList tlist2;
     145
     146  plist2.AddToList(&tlist2);
     147
     148  MReadMarsFile read2("OriginalMC");
     149  read2.AddFile(AnalysisFilename);
     150  read2.DisableAutoScheme();
     151
     152  tlist2.AddToList(&read2);
     153
     154  MWriteRootFile writeOrig(OutFilename->Data(),"UPDATE");
     155  writeOrig.AddContainer("MMcEvtBasic", "OriginalMC", kFALSE);
     156
     157  tlist2.AddToList(&writeOrig);
     158
     159  MEvtLoop evtloop2;
     160  bar.SetWindowName("Copying OriginalMC tree...");
     161  evtloop2.SetProgressBar(&bar);
     162  evtloop2.SetParList(&plist2);
     163 
     164  if (!evtloop2.Eventloop())
     165    return;
     166
     167  tlist2.PrintStatistics();
     168
    180169  return;
    181170}
Note: See TracChangeset for help on using the changeset viewer.