Changeset 5208 for trunk


Ignore:
Timestamp:
10/08/04 10:13:04 (20 years ago)
Author:
rico
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars/mtemp/mifae
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • trunk/MagicSoft/Mars/mtemp/mifae/Changelog

    r5191 r5208  
    1919                                                 -*-*- END OF LINE -*-*-
    2020
     21  2004/10/08  Javier Rico
     22    * programs/makeHillas.cc, programs/makehillas.datacard
     23     - Include possibility to select signal extractor in datacard
     24     - Include possitility to compute pedestals from data themselves,
     25       by selecting word PRUNS 0. For the time being, the old option
     26       of computing them from a pedestal run is still available
     27       
    2128  2004/10/05  Javier Rico
    2229    * programs/srcPos.cc

  • trunk/MagicSoft/Mars/mtemp/mifae/programs/makeHillas.cc

    r5191 r5208  
    2222#include "MExtractor.h"
    2323#include "MExtractFixedWindow.h"
     24#include "MExtractFixedWindowPeakSearch.h"
    2425#include "MExtractSlidingWindow.h"
     26#include "MPedCalcFromLoGain.h"
    2527#include "MExtractSignal.h"
    2628#include "MCalibrationChargeCalc.h"
    … …  
    3032#include "MCerPhotEvt.h"
    3133#include "MPedPhotCam.h"
    32 #include "MCalibrate.h"
     34#include "MCalibrateData.h"
    3335#include "MPedPhotCalc.h"
    3436#include "MHillas.h"
    … …  
    6567void makeHillas();
    6668
    67 // initial and final time slices to be used in signal extraction
    68 const Byte_t hifirst = 1;
    69 const Byte_t hilast  = 14;
    70 const Byte_t lofirst = 3;
    71 const Byte_t lolast  = 14;
    7269
    7370// declaration of variables read from datacards
    … …  
    8481ULong_t  nmaxevents= 999999999;
    8582Short_t  calflag   = 1;
    86 Bool_t  caltimeflag   = kFALSE;
     83Bool_t   caltimeflag   = kFALSE;
    8784Short_t  cleanflag   = 1;
    88 UShort_t  lrings     = 1;
     85UShort_t lrings     = 1;
    8986Float_t  lcore     = 3.0;
    9087Float_t  ltail     = 1.5;
    … …  
    9491Int_t    kalgorithm = 1;
    9592Int_t    nfiles    = 0;
     93Int_t    hifirst = 0;
     94Int_t    hilast  = 14;
     95Int_t    lofirst = 2;
     96Int_t    lolast  = 14;
     97Int_t    wsize   = 6;
     98Int_t    sext    = 0;
    9699
    97100const TString defaultcard="makehillas.datacard";
     101char* chext[3]={"Fixed window","Sliding window","Peak Search"};
    98102/*************************************************************/
    99103static void Usage()
    … …  
    102106  gLog << "Usage is:" << endl;
    103107  gLog << "   makeHillas [-h] [-?] <datacards>" << endl << endl;
    104   gLog << "     <datacards>: datacards file name (dafault input.datacards)" << endl;
     108  gLog << "     <datacards>: datacards file name (dafault makehillas.datacards)" << endl;
    105109  gLog << "     -?/-h: This help" << endl << endl;
    106110}
    … …  
    135139void makeHillas()
    136140{
    137   // Set the general tasks/containers
    138   MExtractFixedWindow    extractor;
    139   extractor.SetRange(hifirst,hilast,lofirst,lolast);
     141  // Set the signal extractor
     142  MExtractor* extractor;
     143  switch(sext)
     144    {
     145    case 0:
     146      extractor = new MExtractFixedWindow();
     147      break;
     148    case 1:
     149      extractor = new MExtractSlidingWindow();
     150      ((MExtractSlidingWindow*)extractor)->SetWindowSize(wsize,wsize);
     151      break;
     152    case 2:
     153      extractor = new MExtractFixedWindowPeakSearch();
     154      ((MExtractFixedWindowPeakSearch*)extractor)->SetWindows(wsize,wsize,4);
     155      break;
     156    default:
     157      extractor = new MExtractFixedWindow();
     158      break;
     159    }
     160
     161  extractor->SetRange(hifirst,hilast,lofirst,lolast);
    140162 
    141   //MExtractSlidingWindow    extractor;
    142   //extractor.SetRange(hifirst,hilast,lofirst,lolast);
    143   //extractor.SetWindowSize(4,4);
    144  
    145   MGeomCamMagic       geomcam;
    146   MGeomApply          geomapl;
    147 
    148   /************************************/
    149   /* FIRST LOOP: PEDESTAL COMPUTATION */
    150   /************************************/
     163
     164  /****************************************************/
     165  /* FIRST LOOP: PEDESTAL COMPUTATION FOR CALIBRATION */
     166  /****************************************************/
    151167
    152168  // If you want to exclude pixels from the beginning, read
    … …  
    155171
    156172  MJPedestal pedcalloop;
    157   pedcalloop.SetInput(&pediter);
    158   pedcalloop.SetExtractor(&extractor);
     173  pedcalloop.SetInput(&pedcaliter);
     174  pedcalloop.SetExtractor(extractor);
    159175  //  pedloop.SetBadPixels(badcam);
    160176  //pedcalloop.SetOutputPath(outpath.Data());
    … …  
    171187
    172188  calloop.SetRelTimeCalibration(caltimeflag);
    173   calloop.SetExtractor(&extractor);
     189  calloop.SetExtractor(extractor);
    174190  calloop.SetInput(&caliter);
    175191  calloop.SetBadPixels(pedcalloop.GetBadPixels());
    … …  
    180196      return;
    181197
    182 
    183   /************************************************************************/
    184   /*                THIRD LOOP: PEDESTAL CALIBRATION INTO PHOTONS         */
    185   /************************************************************************/
    186   // First Compute the pedestals
    187   MJPedestal pedloop;
    188   pedloop.SetInput(&pediter);
    189 
    190   if (!pedloop.Process())
    191     return;
    192 
    193   MParList  plist3;
    194   MTaskList tlist3;
    195   plist3.AddToList(&tlist3);
    196  
    197   // containers
    198   MCerPhotEvt         nphot;
    199   MPedPhotCam         nphotrms;
    200   MExtractedSignalCam sigcam;
    201 
    202   plist3.AddToList(&geomcam);
    203   plist3.AddToList(&pedloop.GetPedestalCam());
    204   plist3.AddToList(&calloop.GetCalibrationCam());
    205   plist3.AddToList(&calloop.GetQECam());
    206   plist3.AddToList(&calloop.GetRelTimeCam());
    207   plist3.AddToList(&calloop.GetBadPixels());
    208   plist3.AddToList(&sigcam);
    209   plist3.AddToList(&nphot);
    210   plist3.AddToList(&nphotrms);
    211 
    212  
    213   MCalibrate::CalibrationMode_t calMode=MCalibrate::kDefault; 
    214   if(calflag==0)
    215     calMode=MCalibrate::kNone;
    216   if(calflag==-1)
    217     calMode=MCalibrate::kDummy;
    218 
    219   //tasks
    220   MReadMarsFile read3("Events");
    221   static_cast<MRead&>(read3).AddFiles(pediter);
    222   read3.DisableAutoScheme();
    223  
    224   MCalibrate      photcalc(calMode);
    225   MPedPhotCalc    photrmscalc;
    226  
    227   tlist3.AddToList(&read3);
    228   tlist3.AddToList(&geomapl);
    229   tlist3.AddToList(&extractor);
    230   tlist3.AddToList(&photcalc);
    231   tlist3.AddToList(&photrmscalc);
    232 
    233   // Create and setup the eventloop
    234   MEvtLoop evtloop3;
    235   evtloop3.SetParList(&plist3);
    236   if (!evtloop3.Eventloop())
    237     return;
    238  
    239   tlist3.PrintStatistics();
    240 
    241   /************************************************************************/
    242   /*                FOURTH LOOP: DATA CALIBRATION INTO PHOTONS            */
    243   /************************************************************************/
    244 
    245   MParList  plist4;
    246   MTaskList tlist4;
    247   plist4.AddToList(&tlist4);
    248  
    249   // containers
    250   MHillas       hillas;
    251   MNewImagePar  newimagepar;
    252   MSrcPosCam    source;
    253   MRawRunHeader runhead;
    254 
    255   MArrivalTimeCam   timecam;
    256 
    257 
    258   // islands
    259   MIslands      isl;
    260   MIslands      isl2;
    261   MIslands      isl3;
    262 
    263   isl.SetName("MIslands"); 
    264   isl2.SetName("MIslands2");
    265   isl3.SetName("MIslands3");
    266 
    267   plist4.AddToList(&timecam);
    268   plist4.AddToList(&isl);
    269  
    270   if (islflag == 2)
    271     plist4.AddToList(&isl2);
    272   if (islflag == 3)
    273     plist4.AddToList(&isl3);
    274  
    275   plist4.AddToList(&geomcam);
    276   plist4.AddToList(&pedloop.GetPedestalCam());
    277   plist4.AddToList(&calloop.GetCalibrationCam());
    278   plist4.AddToList(&calloop.GetQECam());
    279   plist4.AddToList(&calloop.GetRelTimeCam());
    280   plist4.AddToList(&calloop.GetBadPixels());
    281   plist4.AddToList(&nphot);
    282   plist4.AddToList(&nphotrms);
    283   plist4.AddToList(&source);
    284   plist4.AddToList(&hillas);
    285   plist4.AddToList(&newimagepar);
    286   plist4.AddToList(&runhead);
    287 
    288   // cuts
    289   MF cut(filter);
    290  
    291   //tasks
    292   MReadMarsFile read4("Events");
    293   static_cast<MRead&>(read4).AddFiles(datiter);
    294   read4.DisableAutoScheme();
    295  
    296   MImgCleanStd      clean(lcore,ltail);
    297   clean.SetCleanRings(lrings);
    298   MImgCleanStd::CleaningMethod_t cleanMeth= MImgCleanStd::kStandard; 
    299   if(cleanflag==2)
    300     cleanMeth=MImgCleanStd::kDemocratic;
    301   clean.SetMethod(cleanMeth);
    302 
    303   MArrivalTimeCalc2 timecalc;
    304   MIslandsCalc       island;
    305   island.SetOutputName("MIslands");
    306   island.SetAlgorithm(kalgorithm);
    307 
    308   MBadPixelsTreat   interpolatebadpixels;
    309   interpolatebadpixels.SetUseInterpolation();
    310   interpolatebadpixels.SetProcessPedestal();
    311   //  interpolatebadpixels.SetSloppyTreatment();
    312 
    313   MIslandsClean      islclean(lnew);
    314   islclean.SetInputName("MIslands");
    315   islclean.SetMethod(kmethod);
    316      
    317   MIslandsCalc       island2;
    318   island2.SetOutputName("MIslands2"); 
    319   island2.SetAlgorithm(kalgorithm);
    320 
    321   MIslandsCalc       island3;
    322   island3.SetOutputName("MIslands3"); 
    323  
    324  
    325   MHillasCalc       hcalc;
    326   MHillasSrcCalc    csrc1;
    327  
    328   MContinue applycut(&cut);
    329   applycut.SetInverted(kTRUE);
    330   MWriteRootFile write(outname,"RECREATE");
    331  
    332   MHillasDisplay*  disphillas=NULL;
    333 
    334   write.AddContainer("MHillas"        , "Parameters");
    335   write.AddContainer("MHillasSrc"     , "Parameters");
    336   write.AddContainer("MHillasExt"     , "Parameters");
    337   write.AddContainer("MNewImagePar"   , "Parameters");
    338   write.AddContainer("MRawEvtHeader"  , "Parameters");
    339   write.AddContainer("MRawRunHeader"  , "Parameters");
    340   write.AddContainer("MTime"          , "Parameters");
    341   write.AddContainer("MConcentration" , "Parameters");
    342   write.AddContainer("MSrcPosCam"     , "Parameters");
    343   write.AddContainer("MIslands"       , "Parameters");
    344 
    345   if (islflag == 2)
    346     write.AddContainer("MIslands2" , "Parameters");
    347   if (islflag == 3)
    348     write.AddContainer("MIslands3" , "Parameters");
    349 
    350 
    351   if(display)
    352     {
    353       disphillas = new MHillasDisplay(&nphot,&geomcam);
    354       disphillas->SetIslandsName("MIslands");
     198  // Next loops are different if we take pedestals from data or from pedestal files
     199  if(pediter.GetNumRuns()==0)
     200    {
     201      /************************************************************************/
     202      /*                 THIRD (SMALL) LOOP TO GET INITIAl PEDESTALS          */
     203      /************************************************************************/     
     204      MParList  plist3;
     205      MTaskList tlist3;
     206     
     207      // containers
     208      MPedestalCam  pedcamdata;     
     209      MGeomCamMagic geomcam;
     210     
     211      plist3.AddToList(&tlist3);
     212      plist3.AddToList(&geomcam);
     213      plist3.AddToList(&pedcamdata);
     214           
     215      //tasks
     216      MReadMarsFile read3("Events");
     217      static_cast<MRead&>(read3).AddFiles(datiter);
     218      read3.DisableAutoScheme();
     219
     220      MGeomApply      geomapl;
     221
     222      MPedCalcFromLoGain peddatacalc;
     223      peddatacalc.SetPedestalUpdate(kFALSE);
     224      peddatacalc.SetRange(hifirst,hilast,lofirst,lolast);
     225      peddatacalc.SetWindowSize(wsize,wsize);
     226      peddatacalc.SetNumEventsDump(500);
     227      peddatacalc.SetMaxHiGainVar(40);
     228     
     229      tlist3.AddToList(&read3);
     230      tlist3.AddToList(&geomapl);
     231      tlist3.AddToList(&peddatacalc);
     232     
     233      // Create and setup the eventloop
     234      MEvtLoop evtloop3;
     235      evtloop3.SetParList(&plist3);
     236      if (!evtloop3.Eventloop(500))
     237        return;
     238     
     239      tlist3.PrintStatistics();
     240     
     241      /************************************************************************/
     242      /*       FOURTH LOOP: PEDESTAL+DATA CALIBRATION INTO PHOTONS            */
     243      /************************************************************************/
     244     
     245      MParList  plist4;
     246      MTaskList tlist4;
     247      plist4.AddToList(&tlist4);
     248     
     249      // containers
     250      MCerPhotEvt     nphot;
     251      MPedPhotCam     nphotrms;
     252      MHillas         hillas;
     253      MNewImagePar    newimagepar;
     254      MSrcPosCam      source;
     255      MRawRunHeader   runhead;     
     256      MArrivalTimeCam timecam;     
     257      // islands
     258      MIslands      isl;
     259      MIslands      isl2;
     260      MIslands      isl3;
     261     
     262      isl.SetName("MIslands"); 
     263      isl2.SetName("MIslands2");
     264      isl3.SetName("MIslands3");
     265     
     266      plist4.AddToList(&timecam);
     267      plist4.AddToList(&isl);
     268     
    355269      if (islflag == 2)
    356         disphillas->SetIslandsName("MIslands2");
     270        plist4.AddToList(&isl2);
    357271      if (islflag == 3)
    358         disphillas->SetIslandsName("MIslands3");
    359     }     
    360 
    361   tlist4.AddToList(&read4);
    362   tlist4.AddToList(&geomapl);
    363   tlist4.AddToList(&extractor);
    364   tlist4.AddToList(&photcalc);
    365   if(calflag==11 || calflag==21)
    366     tlist4.AddToList(&interpolatebadpixels);
    367   tlist4.AddToList(&clean);
    368   tlist4.AddToList(&timecalc);
    369   tlist4.AddToList(&island);
    370 
    371   if (islflag == 2)
    372     {
    373       tlist4.AddToList(&islclean);
    374       tlist4.AddToList(&island2);
    375     }
    376 
    377   if (islflag == 3)
    378     {
    379       tlist4.AddToList(&islclean);
    380       tlist4.AddToList(&island3);
    381     }
    382  
    383   tlist4.AddToList(&hcalc);
    384   tlist4.AddToList(&csrc1);
    385   if(filter.Length())
    386     tlist4.AddToList(&applycut);
    387   tlist4.AddToList(&write);
    388   if(display)
    389     {
    390       disphillas->SetPSFile();
    391       disphillas->SetPSFileName(psfilename);
    392       if(display==2)
    393         disphillas->SetPause(kFALSE);   
    394       tlist4.AddToList(disphillas);
    395     }
    396 
    397   // Create and setup the eventloop
    398   MEvtLoop datloop;
    399   datloop.SetParList(&plist4);
    400 
    401   cout << "*************************************************************" << endl;
    402   cout << "***   COMPUTING DATA USING EXTRACTED SIGNAL (IN PHOTONS)  ***" << endl;
    403   cout << "*************************************************************" << endl;
    404  
    405   if (!datloop.Eventloop(nmaxevents))
    406     return;
    407 
    408   tlist4.PrintStatistics();   
    409 
     272        plist4.AddToList(&isl3);
     273     
     274      plist4.AddToList(&geomcam);
     275      plist4.AddToList(&pedcamdata);
     276      plist4.AddToList(&calloop.GetCalibrationCam());
     277      plist4.AddToList(&calloop.GetQECam());
     278      plist4.AddToList(&calloop.GetRelTimeCam());
     279      plist4.AddToList(&calloop.GetBadPixels());
     280      plist4.AddToList(&nphot);
     281      plist4.AddToList(&nphotrms);
     282      plist4.AddToList(&source);
     283      plist4.AddToList(&hillas);
     284      plist4.AddToList(&newimagepar);
     285      plist4.AddToList(&runhead);
     286     
     287      // cuts
     288      MF cut(filter);
     289     
     290      //tasks
     291      MReadMarsFile read4("Events");
     292      static_cast<MRead&>(read4).AddFiles(datiter);
     293      read4.DisableAutoScheme();
     294           
     295      peddatacalc.SetPedestalUpdate(kTRUE);
     296
     297      MCalibrateData::CalibrationMode_t calMode=MCalibrateData::kDefault; 
     298      if(calflag==0)
     299        calMode=MCalibrateData::kNone;
     300      if(calflag==-1)
     301        calMode=MCalibrateData::kDummy;
     302      MCalibrateData  photcalc;
     303      photcalc.SetCalibrationMode(calMode);
     304      photcalc.EnablePedestalType(MCalibrateData::kRun);
     305      photcalc.EnablePedestalType(MCalibrateData::kEvent);
     306
     307      MImgCleanStd      clean(lcore,ltail);
     308      clean.SetCleanRings(lrings);
     309      MImgCleanStd::CleaningMethod_t cleanMeth= MImgCleanStd::kStandard; 
     310      if(cleanflag==2)
     311        cleanMeth=MImgCleanStd::kDemocratic;
     312      clean.SetMethod(cleanMeth);
     313     
     314      MArrivalTimeCalc2 timecalc;
     315      MIslandsCalc       island;
     316      island.SetOutputName("MIslands");
     317      island.SetAlgorithm(kalgorithm);
     318     
     319      MBadPixelsTreat   interpolatebadpixels;
     320      interpolatebadpixels.SetUseInterpolation();
     321      interpolatebadpixels.SetProcessPedestal();
     322     
     323      MIslandsClean      islclean(lnew);
     324      islclean.SetInputName("MIslands");
     325      islclean.SetMethod(kmethod);
     326     
     327      MIslandsCalc       island2;
     328      island2.SetOutputName("MIslands2"); 
     329      island2.SetAlgorithm(kalgorithm);
     330     
     331      MIslandsCalc       island3;
     332      island3.SetOutputName("MIslands3"); 
     333     
     334     
     335      MHillasCalc       hcalc;
     336      MHillasSrcCalc    csrc1;
     337     
     338      MContinue applycut(&cut);
     339      applycut.SetInverted(kTRUE);
     340      MWriteRootFile write(outname,"RECREATE");
     341     
     342      MHillasDisplay*  disphillas=NULL;
     343     
     344      write.AddContainer("MHillas"        , "Parameters");
     345      write.AddContainer("MHillasSrc"     , "Parameters");
     346      write.AddContainer("MHillasExt"     , "Parameters");
     347      write.AddContainer("MNewImagePar"   , "Parameters");
     348      write.AddContainer("MRawEvtHeader"  , "Parameters");
     349      write.AddContainer("MRawRunHeader"  , "Parameters");
     350      write.AddContainer("MTime"          , "Parameters");
     351      write.AddContainer("MConcentration" , "Parameters");
     352      write.AddContainer("MSrcPosCam"     , "Parameters");
     353      write.AddContainer("MIslands"       , "Parameters");
     354     
     355      if (islflag == 2)
     356        write.AddContainer("MIslands2" , "Parameters");
     357      if (islflag == 3)
     358        write.AddContainer("MIslands3" , "Parameters");
     359     
     360     
     361      if(display)
     362        {
     363          disphillas = new MHillasDisplay(&nphot,&geomcam);
     364          disphillas->SetIslandsName("MIslands");
     365          if (islflag == 2)
     366            disphillas->SetIslandsName("MIslands2");
     367          if (islflag == 3)
     368            disphillas->SetIslandsName("MIslands3");
     369        }     
     370     
     371      tlist4.AddToList(&read4);
     372      tlist4.AddToList(&geomapl);
     373      tlist4.AddToList(&peddatacalc);
     374      tlist4.AddToList(extractor);
     375      tlist4.AddToList(&timecalc);
     376      tlist4.AddToList(&photcalc);
     377      if(calflag==11 || calflag==21)
     378        tlist4.AddToList(&interpolatebadpixels);
     379      tlist4.AddToList(&clean);
     380      tlist4.AddToList(&island);
     381     
     382      if (islflag == 2)
     383        {
     384          tlist4.AddToList(&islclean);
     385          tlist4.AddToList(&island2);
     386        }
     387     
     388      if (islflag == 3)
     389        {
     390          tlist4.AddToList(&islclean);
     391          tlist4.AddToList(&island3);
     392        }
     393     
     394      tlist4.AddToList(&hcalc);
     395      tlist4.AddToList(&csrc1);
     396      if(filter.Length())
     397        tlist4.AddToList(&applycut);
     398      tlist4.AddToList(&write);
     399      if(display)
     400        {
     401          disphillas->SetPSFile();
     402          disphillas->SetPSFileName(psfilename);
     403          if(display==2)
     404            disphillas->SetPause(kFALSE);       
     405          tlist4.AddToList(disphillas);
     406        }
     407     
     408      // Create and setup the eventloop
     409      MEvtLoop datloop;
     410      datloop.SetParList(&plist4);
     411     
     412      cout << "*************************************************************" << endl;
     413      cout << "***   COMPUTING DATA USING EXTRACTED SIGNAL (IN PHOTONS)  ***" << endl;
     414      cout << "*************************************************************" << endl;
     415     
     416      if (!datloop.Eventloop(nmaxevents))
     417        return;
     418     
     419      tlist4.PrintStatistics();   
     420      delete extractor;
     421    }
     422  else
     423    {     
     424      /************************************************************************/
     425      /*                THIRD LOOP: PEDESTAL CALIBRATION INTO PHOTONS         */
     426      /************************************************************************/
     427     
     428      // First Compute the pedestals
     429      MJPedestal pedloop;
     430      pedloop.SetInput(&pediter);
     431     
     432      if (!pedloop.Process())
     433        return;
     434     
     435      MParList  plist3;
     436      MTaskList tlist3;
     437      plist3.AddToList(&tlist3);
     438     
     439      // containers
     440      MGeomCamMagic       geomcam;
     441      MCerPhotEvt         nphot;
     442      MPedPhotCam         nphotrms;
     443      MExtractedSignalCam sigcam;
     444     
     445      plist3.AddToList(&geomcam);
     446      plist3.AddToList(&pedloop.GetPedestalCam());
     447      plist3.AddToList(&calloop.GetCalibrationCam());
     448      plist3.AddToList(&calloop.GetQECam());
     449      plist3.AddToList(&calloop.GetRelTimeCam());
     450      plist3.AddToList(&calloop.GetBadPixels());
     451      plist3.AddToList(&sigcam);
     452      plist3.AddToList(&nphot);
     453      plist3.AddToList(&nphotrms);
     454     
     455     
     456      MCalibrateData::CalibrationMode_t calMode=MCalibrateData::kDefault; 
     457      if(calflag==0)
     458        calMode=MCalibrateData::kNone;
     459      if(calflag==-1)
     460        calMode=MCalibrateData::kDummy;
     461     
     462      //tasks
     463      MReadMarsFile read3("Events");
     464      static_cast<MRead&>(read3).AddFiles(pediter);
     465      read3.DisableAutoScheme();
     466     
     467      MGeomApply      geomapl;
     468      MPedPhotCalc    photrmscalc;
     469      MCalibrateData  photcalc;
     470      photcalc.SetCalibrationMode(calMode);
     471     
     472      tlist3.AddToList(&read3);
     473      tlist3.AddToList(&geomapl);
     474      tlist3.AddToList(extractor);
     475      tlist3.AddToList(&photcalc);
     476      tlist3.AddToList(&photrmscalc);
     477     
     478      // Create and setup the eventloop
     479      MEvtLoop evtloop3;
     480      evtloop3.SetParList(&plist3);
     481      if (!evtloop3.Eventloop())
     482        return;
     483     
     484      tlist3.PrintStatistics();
     485     
     486      /************************************************************************/
     487      /*                FOURTH LOOP: DATA CALIBRATION INTO PHOTONS            */
     488      /************************************************************************/
     489     
     490      MParList  plist4;
     491      MTaskList tlist4;
     492      plist4.AddToList(&tlist4);
     493     
     494      // containers
     495      MHillas       hillas;
     496      MNewImagePar  newimagepar;
     497      MSrcPosCam    source;
     498      MRawRunHeader runhead;
     499     
     500      MArrivalTimeCam   timecam;
     501     
     502     
     503      // islands
     504      MIslands      isl;
     505      MIslands      isl2;
     506      MIslands      isl3;
     507     
     508      isl.SetName("MIslands"); 
     509      isl2.SetName("MIslands2");
     510      isl3.SetName("MIslands3");
     511     
     512      plist4.AddToList(&timecam);
     513      plist4.AddToList(&isl);
     514     
     515      if (islflag == 2)
     516        plist4.AddToList(&isl2);
     517      if (islflag == 3)
     518        plist4.AddToList(&isl3);
     519     
     520      plist4.AddToList(&geomcam);
     521      plist4.AddToList(&pedloop.GetPedestalCam());
     522      plist4.AddToList(&calloop.GetCalibrationCam());
     523      plist4.AddToList(&calloop.GetQECam());
     524      plist4.AddToList(&calloop.GetRelTimeCam());
     525      plist4.AddToList(&calloop.GetBadPixels());
     526      plist4.AddToList(&nphot);
     527      plist4.AddToList(&nphotrms);
     528      plist4.AddToList(&source);
     529      plist4.AddToList(&hillas);
     530      plist4.AddToList(&newimagepar);
     531      plist4.AddToList(&runhead);
     532     
     533      // cuts
     534      MF cut(filter);
     535     
     536      //tasks
     537      MReadMarsFile read4("Events");
     538      static_cast<MRead&>(read4).AddFiles(datiter);
     539      read4.DisableAutoScheme();
     540     
     541      MImgCleanStd      clean(lcore,ltail);
     542      clean.SetCleanRings(lrings);
     543      MImgCleanStd::CleaningMethod_t cleanMeth= MImgCleanStd::kStandard; 
     544      if(cleanflag==2)
     545        cleanMeth=MImgCleanStd::kDemocratic;
     546      clean.SetMethod(cleanMeth);
     547     
     548      MArrivalTimeCalc2 timecalc;
     549      MIslandsCalc       island;
     550      island.SetOutputName("MIslands");
     551      island.SetAlgorithm(kalgorithm);
     552     
     553      MBadPixelsTreat   interpolatebadpixels;
     554      interpolatebadpixels.SetUseInterpolation();
     555      interpolatebadpixels.SetProcessPedestal();
     556      //  interpolatebadpixels.SetSloppyTreatment();
     557     
     558      MIslandsClean      islclean(lnew);
     559      islclean.SetInputName("MIslands");
     560      islclean.SetMethod(kmethod);
     561     
     562      MIslandsCalc       island2;
     563      island2.SetOutputName("MIslands2"); 
     564      island2.SetAlgorithm(kalgorithm);
     565     
     566      MIslandsCalc       island3;
     567      island3.SetOutputName("MIslands3"); 
     568     
     569     
     570      MHillasCalc       hcalc;
     571      MHillasSrcCalc    csrc1;
     572     
     573      MContinue applycut(&cut);
     574      applycut.SetInverted(kTRUE);
     575      MWriteRootFile write(outname,"RECREATE");
     576     
     577      MHillasDisplay*  disphillas=NULL;
     578     
     579      write.AddContainer("MHillas"        , "Parameters");
     580      write.AddContainer("MHillasSrc"     , "Parameters");
     581      write.AddContainer("MHillasExt"     , "Parameters");
     582      write.AddContainer("MNewImagePar"   , "Parameters");
     583      write.AddContainer("MRawEvtHeader"  , "Parameters");
     584      write.AddContainer("MRawRunHeader"  , "Parameters");
     585      write.AddContainer("MTime"          , "Parameters");
     586      write.AddContainer("MConcentration" , "Parameters");
     587      write.AddContainer("MSrcPosCam"     , "Parameters");
     588      write.AddContainer("MIslands"       , "Parameters");
     589     
     590      if (islflag == 2)
     591        write.AddContainer("MIslands2" , "Parameters");
     592      if (islflag == 3)
     593        write.AddContainer("MIslands3" , "Parameters");
     594     
     595     
     596      if(display)
     597        {
     598          disphillas = new MHillasDisplay(&nphot,&geomcam);
     599          disphillas->SetIslandsName("MIslands");
     600          if (islflag == 2)
     601            disphillas->SetIslandsName("MIslands2");
     602          if (islflag == 3)
     603            disphillas->SetIslandsName("MIslands3");
     604        }     
     605     
     606      tlist4.AddToList(&read4);
     607      tlist4.AddToList(&geomapl);
     608      tlist4.AddToList(extractor);
     609      tlist4.AddToList(&photcalc);
     610      if(calflag==11 || calflag==21)
     611        tlist4.AddToList(&interpolatebadpixels);
     612      tlist4.AddToList(&clean);
     613      tlist4.AddToList(&timecalc);
     614      tlist4.AddToList(&island);
     615     
     616      if (islflag == 2)
     617        {
     618          tlist4.AddToList(&islclean);
     619          tlist4.AddToList(&island2);
     620        }
     621     
     622      if (islflag == 3)
     623        {
     624          tlist4.AddToList(&islclean);
     625          tlist4.AddToList(&island3);
     626        }
     627     
     628      tlist4.AddToList(&hcalc);
     629      tlist4.AddToList(&csrc1);
     630      if(filter.Length())
     631        tlist4.AddToList(&applycut);
     632      tlist4.AddToList(&write);
     633      if(display)
     634        {
     635          disphillas->SetPSFile();
     636          disphillas->SetPSFileName(psfilename);
     637          if(display==2)
     638            disphillas->SetPause(kFALSE);       
     639          tlist4.AddToList(disphillas);
     640        }
     641     
     642      // Create and setup the eventloop
     643      MEvtLoop datloop;
     644      datloop.SetParList(&plist4);
     645     
     646      cout << "*************************************************************" << endl;
     647      cout << "***   COMPUTING DATA USING EXTRACTED SIGNAL (IN PHOTONS)  ***" << endl;
     648      cout << "*************************************************************" << endl;
     649     
     650      if (!datloop.Eventloop(nmaxevents))
     651        return;
     652     
     653      tlist4.PrintStatistics();   
     654      delete extractor;
     655    }
    410656}
    411 //-------------------------------------------------------------------------------
     657  //-------------------------------------------------------------------------------
    412658
    413659Bool_t readDatacards(TString& filename)
    … …  
    450696            cout << "readDataCards Warning: adding pedestal runs to the existing list" << endl;
    451697          ifun >> word;
    452           pediter.AddRuns(word.Data(),idirname.Data());
     698          if(strcmp(word.Data(),"0")!=0)
     699            pediter.AddRuns(word.Data(),idirname.Data());
    453700        }
    454701
    … …  
    530777        }
    531778     
     779      // cleaning level
     780      if(strcmp(word.Data(),"EXTRACTOR")==0)
     781        ifun >> sext >> hifirst >> hilast >> lofirst >> lolast >> wsize;
     782     
    532783      if(strcmp(word.Data(),"ISLFLAG")==0)
    533784        {
    … …  
    549800
    550801  pediter.Reset();
     802  pedcaliter.Reset();
    551803  caliter.Reset();
    552804  datiter.Reset();
    … …  
    557809  cout << "* Datacards read from file " << filename << endl;
    558810  cout << "************************************************" << endl;
    559   cout << "Pedestal file (s): "  << endl;
    560   while(!(pfile=pediter.Next()).IsNull())
     811  cout << "Pedestal file (s) for calibration: "  << endl;
     812  while(!(pfile=pedcaliter.Next()).IsNull())
    561813    cout << pfile << endl;
    562814  cout << "Calibration file (s): "  << endl;
    563815  while(!(pfile=caliter.Next()).IsNull())
    564816    cout << pfile << endl;
     817  if(pediter.GetNumRuns()>0)
     818    {
     819      cout << "Pedestal file (s): "  << endl;
     820      while(!(pfile=pediter.Next()).IsNull())
     821        cout << pfile << endl;
     822    }
     823  else
     824    cout << "Warning: Pedestals for data will be computed from data themselves" << endl;
    565825  cout << "Data file (s): "  << endl;
    566826  while(!(pfile=datiter.Next()).IsNull())
    … …  
    572832  if(display)
    573833    cout << "Generating PS file: " << psfilename << endl;
     834  cout << "Signal Extractor: " << chext[sext] << " with bounds ("<<hifirst<<","<<hilast<<","<<lofirst<<","<<lolast<<"), window size " << wsize << endl;
    574835  cout << "Calibration: ";
    575836  if(calflag==0)   
    576     cout << "Pixel area proportional intercalibration" << endl;
     837    cout << "Pixel area proportional intercalibration (kNone)" << endl;
    577838  else if(calflag==-1)   
    578839    cout << "No calibration whatsoever" << endl;
    … …  
    582843    cout << "Default calibration + bad pixels interpolation" << endl;
    583844  cout << "Cleaning level: ("<<lcore<<","<<ltail<<") - " << lrings << "ring" << endl;
    584   cout << "Cleaning methode: "<< cleanflag << endl;
     845  cout << "Cleaning method: "<< cleanflag << endl;
    585846  if (islflag == 1 || islflag == 2)
    586847    cout << "Island calcultation..." << "using algorithm #" << kalgorithm <<endl;
    … …  
    591852  cout << "***********" << endl << endl;
    592853 
    593   if(!pediter.GetNumEntries())
    594     {
    595       cout << "No pedestal file name specified" << endl;
     854  if(!pedcaliter.GetNumEntries())
     855    {
     856      cout << "No pedestal file for calibration specified" << endl;
    596857      return kFALSE;
    597858    }

  • trunk/MagicSoft/Mars/mtemp/mifae/programs/makehillas.datacard

    r4428 r5208  
    1111
    1212// Maximun number of (data) events to be processed)
    13 NEVENTS 9999999
     13NEVENTS 99999999
    1414
    1515// data file directory
    16 IDIR /local_disk/jrico/rootdata/Crab20040215/
     16IDIR /mnt/magic/Data/rootdata/CrabNebula/Period021/2004_09_22/
    1717
    1818// Pedestal (PRUNS), calibration (CRUNS) and data runs (DRUNS), e.g 1500-23444,25444,25455-26544
    19 PRUNS 16743
    20 CRUNS 16744
    21 DRUNS 16745
     19// if PRUNS 0, take pedestals from low gains
     20PCRUNS 39313
     21CRUNS  39315
     22PRUNS 0
     23DRUNS 39261,39262,39264,39265,39269,39271,39273,39274,39276,39277,39279,39280,39282,39283,39285
     24//,39286,39288,39290,39292-39294,39296-39299,39301,39302,39304-39306,39308,39309,39311,39318,39319,39321,39322,39324,39325,39327,39328,39330,39331,39333,39334,39336-39338,39340-39342,39344,39348-39350,39352-39355,39357,39358,39360-39362,39364-39366,39368-39370,39372,39373
     25// 37478-37480,37482-37484,37486,37487,37489-37491,37493-37495,37497-37499,37501-37502,37504-37505,37507-37508,37510-37512,37514,37516-37519,37521,37522,37524,37525,37527-37529,37531,37532,37534,37535,37537-37539,37541-37544,37546-37551,37553-37555
    2226
    2327// output file name
    2428// OUTFILE ~/magic/mars/mars/hillasCrab/crab20040215OnA.root
    25 OUTFILE ./prueba.root
     29OUTFILE /local_disk/jrico/prueba2.root
    2630
    2731// Selection cut.
    … …  
    4044PSFILENAME makehillas.ps
    4145
     46// Signal extractor, higain first, higainlast, logainfirst, logainlast, window size
     47// signal extractors: 0 fixed window, 1 sliding window, 2 peak search
     48EXTRACTOR 2 0 14 2 14 6
     49
    4250// calibration flag:
    4351// -1: kDummy
    … …  
    4755// 11: kDefault(F factor) + bad pixel interpolation
    4856// 21: kDemocratic + bad pixel interpolation
    49 CALFLAG 0
     57CALFLAG 1
    5058
    5159// calibration time:
    5260//   0: kFALSE (no time calibration)
    5361//   1: kTRUE  (time calibration)
    54 //CALTIME 1
     62CALTIME 1
    5563
    5664
    5765// Cleaning level (tail cut, boundary cut, number of rings, cleaning method - 1=standard, 2=democratic)
    58 CLEANLEVEL 3.0 1.5 1 1
     66CLEANLEVEL 2.5 2.0 1 1
    5967
    6068//Island calculations
    61 // 0  nothing about islands    1:normal algorithm
     69// 0  same as 1                1:normal algorithm
    6270// 1  islands w/o cleaning     2:alternative algorithm
    6371// 2  islands  w  cleaning
Note: See TracChangeset for help on using the changeset viewer.